Skip to content
Snippets Groups Projects
Commit 771be901 authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Framework copied from acts examples

parents
Branches
No related tags found
No related merge requests found
Showing
with 1259 additions and 0 deletions
add_subdirectory(Digitization)
add_subdirectory(Fatras)
add_subdirectory(Geant4)
add_subdirectory(Geant4DD4hep)
add_subdirectory(Generators)
add_subdirectory(GeneratorsPythia8)
add_subdirectory(MaterialMapping)
add_subdirectory(Printers)
add_subdirectory(Propagation)
add_subdirectory(TruthTracking)
add_subdirectory(Vertexing)
add_subdirectory(Fitting)
add_subdirectory(TrackFinding)
find_package(Acts REQUIRED COMPONENTS DigitizationPlugin IdentificationPlugin)
add_library(
ActsExamplesDigitization SHARED
src/DigitizationAlgorithm.cpp
src/HitSmearing.cpp)
target_include_directories(
ActsExamplesDigitization
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>)
target_link_libraries(
ActsExamplesDigitization
PRIVATE
ActsCore ActsDigitizationPlugin ActsIdentificationPlugin
ActsExamplesFramework
Boost::program_options)
install(
TARGETS ActsExamplesDigitization
LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR})
// This file is part of the Acts project.
//
// Copyright (C) 2017 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <memory>
#include <string>
#include <unordered_map>
#include "ACTFW/Framework/BareAlgorithm.hpp"
#include "ACTFW/Framework/RandomNumbers.hpp"
#include "Acts/Geometry/GeometryID.hpp"
namespace Acts {
class DigitizationModule;
class IdentifiedDetectorElement;
class PlanarModuleStepper;
class Surface;
class TrackingGeometry;
} // namespace Acts
namespace FW {
/// Create planar clusters from simulation hits.
class DigitizationAlgorithm final : public BareAlgorithm {
public:
struct Config {
/// Input collection of simulated hits.
std::string inputSimulatedHits;
/// Output collection of clusters.
std::string outputClusters;
/// Tracking geometry required to access global-to-local transforms.
std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry;
/// Module stepper for geometric clustering.
std::shared_ptr<const Acts::PlanarModuleStepper> planarModuleStepper;
/// Random numbers tool.
std::shared_ptr<const RandomNumbers> randomNumbers;
};
/// Construct the digitization algorithm.
///
/// @param cfg is the algorithm configuration
/// @param lvl is the logging level
DigitizationAlgorithm(Config cfg, Acts::Logging::Level lvl);
/// Build clusters from input simulation hits.
///
/// @param txt is the algorithm context with event information
/// @return a process code indication success or failure
ProcessCode execute(const AlgorithmContext& ctx) const final override;
private:
struct Digitizable {
const Acts::Surface* surface = nullptr;
const Acts::IdentifiedDetectorElement* detectorElement = nullptr;
const Acts::DigitizationModule* digitizer = nullptr;
};
Config m_cfg;
/// Lookup container for all digitizable surfaces
std::unordered_map<Acts::GeometryID, Digitizable> m_digitizables;
};
} // namespace FW
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <string>
#include <unordered_map>
#include "ACTFW/Framework/BareAlgorithm.hpp"
#include "ACTFW/Framework/RandomNumbers.hpp"
#include "Acts/Geometry/GeometryID.hpp"
namespace Acts {
class Surface;
class TrackingGeometry;
} // namespace Acts
namespace FW {
/// Create fittable measurements using truth smearing.
///
/// The truth information is smeared in the local measurement frame using
/// Gaussian noise to generate a fittable measurement, i.e. a source link.
class HitSmearing final : public BareAlgorithm {
public:
struct Config {
/// Input collection of simulated hits.
std::string inputSimulatedHits;
/// Output collection for source links with smeared measurements.
std::string outputSourceLinks;
/// Width of the Gaussian smearing, i.e. resolution; must be positive.
double sigmaLoc0 = -1;
double sigmaLoc1 = -1;
/// Tracking geometry required to access global-to-local transforms.
std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry;
/// Random numbers tool.
std::shared_ptr<const RandomNumbers> randomNumbers = nullptr;
};
HitSmearing(const Config& cfg, Acts::Logging::Level lvl);
ProcessCode execute(const AlgorithmContext& ctx) const final override;
private:
Config m_cfg;
/// Lookup container for hit surfaces that generate smeared hits
std::unordered_map<Acts::GeometryID, const Acts::Surface*> m_surfaces;
};
} // namespace FW
// This file is part of the Acts project.
//
// Copyright (C) 2017 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "ACTFW/Digitization/DigitizationAlgorithm.hpp"
#include <iostream>
#include <stdexcept>
#include "ACTFW/EventData/GeometryContainers.hpp"
#include "ACTFW/EventData/SimHit.hpp"
#include "ACTFW/EventData/SimParticle.hpp"
#include "ACTFW/EventData/SimVertex.hpp"
#include "ACTFW/Framework/WhiteBoard.hpp"
#include "Acts/EventData/TrackParameters.hpp"
#include "Acts/Geometry/DetectorElementBase.hpp"
#include "Acts/Geometry/GeometryID.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/Plugins/Digitization/DigitizationModule.hpp"
#include "Acts/Plugins/Digitization/PlanarModuleCluster.hpp"
#include "Acts/Plugins/Digitization/PlanarModuleStepper.hpp"
#include "Acts/Plugins/Digitization/Segmentation.hpp"
#include "Acts/Plugins/Identification/IdentifiedDetectorElement.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/ParameterDefinitions.hpp"
#include "Acts/Utilities/Units.hpp"
FW::DigitizationAlgorithm::DigitizationAlgorithm(
FW::DigitizationAlgorithm::Config cfg, Acts::Logging::Level lvl)
: FW::BareAlgorithm("DigitizationAlgorithm", lvl), m_cfg(std::move(cfg)) {
if (m_cfg.inputSimulatedHits.empty()) {
throw std::invalid_argument("Missing input hits collection");
}
if (m_cfg.outputClusters.empty()) {
throw std::invalid_argument("Missing output clusters collection");
}
if (not m_cfg.trackingGeometry) {
throw std::invalid_argument("Missing tracking geometry");
}
if (!m_cfg.planarModuleStepper) {
throw std::invalid_argument("Missing planar module stepper");
}
if (!m_cfg.randomNumbers) {
throw std::invalid_argument("Missing random numbers tool");
}
// fill the digitizables map to allow lookup by geometry id only
m_cfg.trackingGeometry->visitSurfaces([this](const Acts::Surface* surface) {
Digitizable dg;
// require a valid surface
dg.surface = surface;
if (not dg.surface) {
return;
}
// require an associated detector element
dg.detectorElement = dynamic_cast<const Acts::IdentifiedDetectorElement*>(
dg.surface->associatedDetectorElement());
if (not dg.detectorElement) {
return;
}
// require an associated digitization module
dg.digitizer = dg.detectorElement->digitizationModule().get();
if (not dg.digitizer) {
return;
}
// record all valid surfaces
this->m_digitizables.insert_or_assign(surface->geoID(), dg);
});
}
FW::ProcessCode FW::DigitizationAlgorithm::execute(
const AlgorithmContext& ctx) const {
// Prepare the input and output collections
const auto& hits =
ctx.eventStore.get<SimHitContainer>(m_cfg.inputSimulatedHits);
FW::GeometryIdMultimap<Acts::PlanarModuleCluster> clusters;
for (auto&& [moduleGeoId, moduleHits] : groupByModule(hits)) {
// can only digitize hits on digitizable surfaces
const auto it = m_digitizables.find(moduleGeoId);
if (it == m_digitizables.end()) {
continue;
}
const auto& dg = it->second;
// local intersection / direction
const auto invTransfrom = dg.surface->transform(ctx.geoContext).inverse();
// use iterators manually so we can retrieve the hit index in the container
for (auto ih = moduleHits.begin(); ih != moduleHits.end(); ++ih) {
const auto& hit = *ih;
const auto idx = hits.index_of(ih);
Acts::Vector2D localIntersect = (invTransfrom * hit.position()).head<2>();
Acts::Vector3D localDirection =
invTransfrom.linear() * hit.unitDirection();
// compute digitization steps
const auto thickness = dg.detectorElement->thickness();
const auto lorentzAngle = dg.digitizer->lorentzAngle();
auto lorentzShift = thickness * std::tan(lorentzAngle);
lorentzShift *= -(dg.digitizer->readoutDirection());
// now calculate the steps through the silicon
std::vector<Acts::DigitizationStep> dSteps =
m_cfg.planarModuleStepper->cellSteps(ctx.geoContext, *dg.digitizer,
localIntersect, localDirection);
// everything under threshold or edge effects
if (!dSteps.size()) {
ACTS_VERBOSE("No steps returned from stepper.");
continue;
}
// lets create a cluster - centroid method
double localX = 0.;
double localY = 0.;
double totalPath = 0.;
// the cells to be used
std::vector<Acts::DigitizationCell> usedCells;
usedCells.reserve(dSteps.size());
// loop over the steps
for (auto dStep : dSteps) {
// @todo implement smearing
localX += dStep.stepLength * dStep.stepCellCenter.x();
localY += dStep.stepLength * dStep.stepCellCenter.y();
totalPath += dStep.stepLength;
usedCells.push_back(Acts::DigitizationCell(dStep.stepCell.channel0,
dStep.stepCell.channel1,
dStep.stepLength));
}
// divide by the total path
localX /= totalPath;
localX += lorentzShift;
localY /= totalPath;
// get the segmentation & find the corresponding cell id
const Acts::Segmentation& segmentation = dg.digitizer->segmentation();
auto binUtility = segmentation.binUtility();
Acts::Vector2D localPosition(localX, localY);
// @todo remove unneccesary conversion
// size_t bin0 = binUtility.bin(localPosition, 0);
// size_t bin1 = binUtility.bin(localPosition, 1);
// size_t binSerialized = binUtility.serialize({{bin0, bin1, 0}});
// the covariance is currently set to 0.
Acts::ActsSymMatrixD<3> cov;
cov << 0.05, 0., 0., 0., 0.05, 0., 0., 0.,
900. * Acts::UnitConstants::ps * Acts::UnitConstants::ps;
// create the planar cluster
Acts::PlanarModuleCluster pCluster(
dg.surface->getSharedPtr(), Identifier(identifier_type(idx), {idx}),
std::move(cov), localX, localY, hit.time(), std::move(usedCells));
// insert into the cluster container. since the input data is already
// sorted by geoId, we should always be able to add at the end.
clusters.emplace_hint(clusters.end(), hit.geometryId(),
std::move(pCluster));
}
}
ACTS_DEBUG("digitized " << hits.size() << " hits into " << clusters.size()
<< " clusters");
// write the clusters to the EventStore
ctx.eventStore.add(m_cfg.outputClusters, std::move(clusters));
return FW::ProcessCode::SUCCESS;
}
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "ACTFW/Digitization/HitSmearing.hpp"
#include "ACTFW/EventData/GeometryContainers.hpp"
#include "ACTFW/EventData/SimHit.hpp"
#include "ACTFW/EventData/SimSourceLink.hpp"
#include "ACTFW/Framework/WhiteBoard.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/Utilities/Definitions.hpp"
FW::HitSmearing::HitSmearing(const Config& cfg, Acts::Logging::Level lvl)
: BareAlgorithm("HitSmearing", lvl), m_cfg(cfg) {
if (m_cfg.inputSimulatedHits.empty()) {
throw std::invalid_argument("Missing input simulated hits collection");
}
if (m_cfg.outputSourceLinks.empty()) {
throw std::invalid_argument("Missing output source links collection");
}
if ((m_cfg.sigmaLoc0 < 0) or (m_cfg.sigmaLoc1 < 0)) {
throw std::invalid_argument("Invalid resolution setting");
}
if (not m_cfg.trackingGeometry) {
throw std::invalid_argument("Missing tracking geometry");
}
if (!m_cfg.randomNumbers) {
throw std::invalid_argument("Missing random numbers tool");
}
// fill the surface map to allow lookup by geometry id only
m_cfg.trackingGeometry->visitSurfaces([this](const Acts::Surface* surface) {
// for now we just require a valid surface
if (not surface) {
return;
}
this->m_surfaces.insert_or_assign(surface->geoID(), surface);
});
}
FW::ProcessCode FW::HitSmearing::execute(const AlgorithmContext& ctx) const {
// setup input and output containers
const auto& hits =
ctx.eventStore.get<SimHitContainer>(m_cfg.inputSimulatedHits);
SimSourceLinkContainer sourceLinks;
sourceLinks.reserve(hits.size());
// setup random number generator
auto rng = m_cfg.randomNumbers->spawnGenerator(ctx);
std::normal_distribution<double> stdNormal(0.0, 1.0);
// setup local covariance
// TODO add support for per volume/layer/module settings
Acts::BoundMatrix cov = Acts::BoundMatrix::Zero();
cov(Acts::eLOC_0, Acts::eLOC_0) = m_cfg.sigmaLoc0 * m_cfg.sigmaLoc0;
cov(Acts::eLOC_1, Acts::eLOC_1) = m_cfg.sigmaLoc1 * m_cfg.sigmaLoc1;
for (auto&& [moduleGeoId, moduleHits] : groupByModule(hits)) {
// check if we should create hits for this surface
const auto is = m_surfaces.find(moduleGeoId);
if (is == m_surfaces.end()) {
continue;
}
// smear all truth hits for this module
const Acts::Surface* surface = is->second;
for (const auto& hit : moduleHits) {
// transform global position into local coordinates
Acts::Vector2D pos(0, 0);
surface->globalToLocal(ctx.geoContext, hit.position(),
hit.unitDirection(), pos);
// smear truth to create local measurement
Acts::BoundVector loc = Acts::BoundVector::Zero();
loc[Acts::eLOC_0] = pos[0] + m_cfg.sigmaLoc0 * stdNormal(rng);
loc[Acts::eLOC_1] = pos[1] + m_cfg.sigmaLoc1 * stdNormal(rng);
// create source link at the end of the container
auto it = sourceLinks.emplace_hint(sourceLinks.end(), *surface, hit, 2,
loc, cov);
// ensure hits and links share the same order to prevent ugly surprises
if (std::next(it) != sourceLinks.end()) {
ACTS_FATAL("The hit ordering broke. Run for your life.");
return ProcessCode::ABORT;
}
}
}
ctx.eventStore.add(m_cfg.outputSourceLinks, std::move(sourceLinks));
return ProcessCode::SUCCESS;
}
add_library(
ActsExamplesFatras SHARED
src/FatrasOptions.cpp)
target_include_directories(
ActsExamplesFatras
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>)
target_link_libraries(
ActsExamplesFatras
PUBLIC ActsCore ActsFatras ActsExamplesFramework Boost::program_options)
install(
TARGETS ActsExamplesFatras
LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR})
// This file is part of the Acts project.
//
// Copyright (C) 2017 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <memory>
#include <string>
#include "ACTFW/EventData/SimHit.hpp"
#include "ACTFW/EventData/SimParticle.hpp"
#include "ACTFW/Framework/BareAlgorithm.hpp"
#include "ACTFW/Framework/RandomNumbers.hpp"
#include "ACTFW/Framework/WhiteBoard.hpp"
namespace FW {
/// Fast track simulation using the Acts propagation and navigation.
///
/// @tparam simulator_t the Fatras simulation kernel type
template <typename simulator_t>
class FatrasAlgorithm final : public BareAlgorithm {
public:
struct Config {
/// The particles input collection.
std::string inputParticles;
/// The simulated particles initial state collection.
std::string outputParticlesInitial;
/// The simulated particles final state collection.
std::string outputParticlesFinal;
/// The simulated hits output collection.
std::string outputHits;
/// The simulator kernel.
simulator_t simulator;
/// Random number service.
std::shared_ptr<const RandomNumbers> randomNumbers;
/// Construct the algorithm config with the simulator kernel.
Config(simulator_t&& simulator_) : simulator(std::move(simulator_)) {}
};
/// Construct the algorithm from a config.
///
/// @param cfg is the configuration struct
/// @param lvl is the logging level
FatrasAlgorithm(Config cfg, Acts::Logging::Level lvl)
: FW::BareAlgorithm("FatrasAlgorithm", lvl), m_cfg(std::move(cfg)) {
ACTS_DEBUG("hits on sensitive surfaces: "
<< m_cfg.simulator.charged.selectHitSurface.sensitive);
ACTS_DEBUG("hits on material surfaces: "
<< m_cfg.simulator.charged.selectHitSurface.material);
ACTS_DEBUG("hits on passive surfaces: "
<< m_cfg.simulator.charged.selectHitSurface.passive);
}
/// Run the simulation for a single event.
///
/// @param ctx the algorithm context containing all event information
FW::ProcessCode execute(const AlgorithmContext& ctx) const final override {
// read input containers
const auto& inputParticles =
ctx.eventStore.get<SimParticleContainer>(m_cfg.inputParticles);
// prepare output containers
SimParticleContainer::sequence_type particlesInitialUnordered;
SimParticleContainer::sequence_type particlesFinalUnordered;
SimHitContainer::sequence_type hitsUnordered;
// reserve appropriate resources
constexpr auto meanHitsPerParticle = 16u;
particlesInitialUnordered.reserve(inputParticles.size());
particlesFinalUnordered.reserve(inputParticles.size());
hitsUnordered.reserve(meanHitsPerParticle * inputParticles.size());
// run the simulation w/ a local random generator
auto rng = m_cfg.randomNumbers->spawnGenerator(ctx);
auto ret = m_cfg.simulator.simulate(
ctx.geoContext, ctx.magFieldContext, rng, inputParticles,
particlesInitialUnordered, particlesFinalUnordered, hitsUnordered);
// fatal error leads to panic
if (not ret.ok()) {
ACTS_FATAL("event " << ctx.eventNumber << " simulation failed with error "
<< ret.error());
return ProcessCode::ABORT;
}
// failed particles are just logged. assumes that failed particles are due
// to edge-cases representing a tiny fraction of the event; not due to a
// fundamental issue.
for (const auto& failed : ret.value()) {
ACTS_ERROR("event " << ctx.eventNumber << " particle " << failed.particle
<< " failed to simulate with error " << failed.error
<< ": " << failed.error.message());
}
// TODO is there a point where too many failed particles or failed particles
// of a particular type (e.g. from hard interaction or any primary
// particle) should also lead to a panic?
ACTS_DEBUG(inputParticles.size() << " input particles");
ACTS_DEBUG(particlesInitialUnordered.size()
<< " simulated particles (initial state)");
ACTS_DEBUG(particlesFinalUnordered.size()
<< " simulated particles (final state)");
ACTS_DEBUG(hitsUnordered.size() << " hits");
// restore ordering for output containers
SimParticleContainer particlesInitial;
SimParticleContainer particlesFinal;
SimHitContainer hits;
particlesInitial.adopt_sequence(std::move(particlesInitialUnordered));
particlesFinal.adopt_sequence(std::move(particlesFinalUnordered));
hits.adopt_sequence(std::move(hitsUnordered));
// store ordered output containers
ctx.eventStore.add(m_cfg.outputParticlesInitial,
std::move(particlesInitial));
ctx.eventStore.add(m_cfg.outputParticlesFinal, std::move(particlesFinal));
ctx.eventStore.add(m_cfg.outputHits, std::move(hits));
return FW::ProcessCode::SUCCESS;
}
private:
Config m_cfg;
};
} // namespace FW
// This file is part of the Acts project.
//
// Copyright (C) 2018-2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <boost/program_options.hpp>
#include <utility>
#include "ACTFW/Fatras/FatrasAlgorithm.hpp"
#include "ACTFW/Utilities/OptionsFwd.hpp"
#include "Acts/Utilities/Logger.hpp"
#include "Acts/Utilities/Units.hpp"
#include "ActsFatras/Kernel/Process.hpp"
#include "ActsFatras/Physics/StandardPhysicsLists.hpp"
namespace FW {
namespace Options {
/// Add Fatras options.
///
/// @param desc The options description to add options to
void addFatrasOptions(Description& desc);
/// Read Fatras options to create the algorithm config.
///
/// @tparam simulator_t type of the simulation kernel
/// @param vars the variables to read from
/// @param simulator the simulation kernel
template <typename simulator_t>
typename FatrasAlgorithm<simulator_t>::Config readFatrasConfig(
const Variables& variables, simulator_t&& simulator) {
using namespace Acts::UnitLiterals;
using Config = typename FatrasAlgorithm<simulator_t>::Config;
using PMin = ActsFatras::Min<ActsFatras::Casts::P>;
ACTS_LOCAL_LOGGER(
Acts::getDefaultLogger("FatrasOptions", Acts::Logging::INFO))
Config cfg(std::forward<simulator_t>(simulator));
// simulation particle cuts both for input and physics lists output
const auto pmin = variables["fatras-pmin-gev"].as<double>() * 1_GeV;
cfg.simulator.selectCharged.template get<PMin>().valMin = pmin;
cfg.simulator.selectNeutral.template get<PMin>().valMin = pmin;
cfg.simulator.charged.physics =
ActsFatras::makeChargedElectroMagneticPhysicsList(pmin);
// all physics process are enabled by default
if (not variables["fatras-em-scattering"].as<bool>()) {
cfg.simulator.charged.physics
.template disable<ActsFatras::detail::StandardScattering>();
}
if (not variables["fatras-em-ionisation"].as<bool>()) {
cfg.simulator.charged.physics
.template disable<ActsFatras::detail::StandardBetheBloch>();
}
if (not variables["fatras-em-radiation"].as<bool>()) {
cfg.simulator.charged.physics
.template disable<ActsFatras::detail::StandardBetheHeitler>();
}
// select hit surfaces for charged particles
const std::string hits = variables["fatras-hits"].as<std::string>();
if (hits == "sensitive") {
ACTS_DEBUG("Configure hits on sensitive surfaces");
cfg.simulator.charged.selectHitSurface.sensitive = true;
cfg.simulator.charged.selectHitSurface.material = false;
cfg.simulator.charged.selectHitSurface.passive = false;
} else if (hits == "material") {
ACTS_DEBUG("Configure hits on material surfaces");
cfg.simulator.charged.selectHitSurface.sensitive = false;
cfg.simulator.charged.selectHitSurface.material = true;
cfg.simulator.charged.selectHitSurface.passive = false;
} else if (hits == "all") {
ACTS_DEBUG("Configure hits on all surfaces");
cfg.simulator.charged.selectHitSurface.sensitive = false;
cfg.simulator.charged.selectHitSurface.material = false;
// this includes sensitive and material surfaces
cfg.simulator.charged.selectHitSurface.passive = true;
} else {
ACTS_WARNING("Invalid hits configuration '" << hits << "'");
// none or unknown type -> record nothing
cfg.simulator.charged.selectHitSurface.sensitive = false;
cfg.simulator.charged.selectHitSurface.material = false;
cfg.simulator.charged.selectHitSurface.passive = false;
}
return cfg;
}
} // namespace Options
} // namespace FW
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "ACTFW/Fatras/FatrasOptions.hpp"
#include <string>
void FW::Options::addFatrasOptions(FW::Options::Description& desc) {
using boost::program_options::bool_switch;
using boost::program_options::value;
auto opt = desc.add_options();
opt("fatras-pmin-gev", value<double>()->default_value(0.5),
"Minimum momentum for simulated particles in GeV");
opt("fatras-em-scattering", value<bool>()->default_value(true),
"Simulate multiple scattering of charged particles");
opt("fatras-em-ionisation", value<bool>()->default_value(true),
"Simulate ionisiation/excitation energy loss of charged particles");
opt("fatras-em-radiation", value<bool>()->default_value(true),
"Simulate radiative energy loss of charged particles");
opt("fatras-hits",
value<std::string>()
->value_name("none|sensitive|material|all")
->default_value("sensitive"),
"Which surfaces should record charged particle hits");
}
add_library(
ActsExamplesFitting SHARED
src/FittingAlgorithm.cpp
#src/FittingAlgorithmFitterFunction.cpp
)
target_include_directories(
ActsExamplesFitting
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>)
target_link_libraries(
ActsExamplesFitting
PUBLIC
ActsCore
ActsExamplesFramework ActsExamplesMagneticField
Boost::program_options)
install(
TARGETS ActsExamplesFitting
LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR})
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <functional>
#include <memory>
#include <vector>
#include "ACTFW/EventData/SimSourceLink.hpp"
#include "ACTFW/EventData/Track.hpp"
#include "ACTFW/Framework/BareAlgorithm.hpp"
#include "ACTFW/Plugins/BField/BFieldOptions.hpp"
#include "Acts/Fitter/KalmanFitter.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
namespace FW {
class FittingAlgorithm final : public BareAlgorithm {
public:
using FitterResult = Acts::Result<Acts::KalmanFitterResult<SimSourceLink>>;
/// Fit function that takes input measurements, initial trackstate and fitter
/// options and returns some fit-specific result.
using FitterFunction = std::function<FitterResult(
const std::vector<SimSourceLink>&, const TrackParameters&,
const Acts::KalmanFitterOptions<Acts::VoidOutlierFinder>&)>;
/// Create the fitter function implementation.
///
/// The magnetic field is intentionally given by-value since the variant
/// contains shared_ptr anyways.
static FitterFunction makeFitterFunction(
std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry,
Options::BFieldVariant magneticField, Acts::Logging::Level lvl);
struct Config {
/// Input source links collection.
std::string inputSourceLinks;
/// Input proto tracks collection, i.e. groups of hit indices.
std::string inputProtoTracks;
/// Input initial track parameter estimates for for each proto track.
std::string inputInitialTrackParameters;
/// Output fitted trajectories collection.
std::string outputTrajectories;
/// Type erased fitter function.
FitterFunction fit;
};
/// Constructor of the fitting algorithm
///
/// @param cfg is the config struct to configure the algorihtm
/// @param level is the logging level
FittingAlgorithm(Config cfg, Acts::Logging::Level lvl);
/// Framework execute method of the fitting algorithm
///
/// @param ctx is the algorithm context that holds event-wise information
/// @return a process code to steer the algporithm flow
FW::ProcessCode execute(const FW::AlgorithmContext& ctx) const final override;
private:
Config m_cfg;
};
} // namespace FW
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "ACTFW/Fitting/FittingAlgorithm.hpp"
#include <stdexcept>
#include "ACTFW/EventData/ProtoTrack.hpp"
#include "ACTFW/EventData/Track.hpp"
#include "ACTFW/Framework/WhiteBoard.hpp"
#include "Acts/Surfaces/PerigeeSurface.hpp"
FW::FittingAlgorithm::FittingAlgorithm(Config cfg, Acts::Logging::Level level)
: FW::BareAlgorithm("FittingAlgorithm", level), m_cfg(std::move(cfg)) {
if (m_cfg.inputSourceLinks.empty()) {
throw std::invalid_argument("Missing input source links collection");
}
if (m_cfg.inputProtoTracks.empty()) {
throw std::invalid_argument("Missing input proto tracks collection");
}
if (m_cfg.inputInitialTrackParameters.empty()) {
throw std::invalid_argument(
"Missing input initial track parameters collection");
}
if (m_cfg.outputTrajectories.empty()) {
throw std::invalid_argument("Missing output trajectories collection");
}
}
FW::ProcessCode FW::FittingAlgorithm::execute(
const FW::AlgorithmContext& ctx) const {
// Read input data
const auto sourceLinks =
ctx.eventStore.get<SimSourceLinkContainer>(m_cfg.inputSourceLinks);
const auto protoTracks =
ctx.eventStore.get<ProtoTrackContainer>(m_cfg.inputProtoTracks);
const auto initialParameters = ctx.eventStore.get<TrackParametersContainer>(
m_cfg.inputInitialTrackParameters);
// Consistency cross checks
if (protoTracks.size() != initialParameters.size()) {
ACTS_FATAL("Inconsistent number of proto tracks and parameters");
return ProcessCode::ABORT;
}
// Prepare the output data with MultiTrajectory
TrajectoryContainer trajectories;
trajectories.reserve(protoTracks.size());
// Construct a perigee surface as the target surface
auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
Acts::Vector3D{0., 0., 0.});
// Perform the fit for each input track
std::vector<SimSourceLink> trackSourceLinks;
for (std::size_t itrack = 0; itrack < protoTracks.size(); ++itrack) {
// The list of hits and the initial start parameters
const auto& protoTrack = protoTracks[itrack];
const auto& initialParams = initialParameters[itrack];
// We can have empty tracks which must give empty fit results
if (protoTrack.empty()) {
trajectories.push_back(SimMultiTrajectory());
ACTS_WARNING("Empty track " << itrack << " found.");
continue;
}
// Clear & reserve the right size
trackSourceLinks.clear();
trackSourceLinks.reserve(protoTrack.size());
// Fill the source links via their indices from the container
for (auto hitIndex : protoTrack) {
auto sourceLink = sourceLinks.nth(hitIndex);
if (sourceLink == sourceLinks.end()) {
ACTS_FATAL("Proto track " << itrack << " contains invalid hit index"
<< hitIndex);
return ProcessCode::ABORT;
}
trackSourceLinks.push_back(*sourceLink);
}
// Set the KalmanFitter options
Acts::KalmanFitterOptions<Acts::VoidOutlierFinder> kfOptions(
ctx.geoContext, ctx.magFieldContext, ctx.calibContext,
Acts::VoidOutlierFinder(), &(*pSurface));
ACTS_DEBUG("Invoke fitter");
auto result = m_cfg.fit(trackSourceLinks, initialParams, kfOptions);
if (result.ok()) {
// Get the fit output object
const auto& fitOutput = result.value();
// The track entry indices container. One element here.
std::vector<size_t> trackTips;
trackTips.reserve(1);
trackTips.emplace_back(fitOutput.trackTip);
// The fitted parameters container. One element (at most) here.
IndexedParams indexedParams;
if (fitOutput.fittedParameters) {
const auto& params = fitOutput.fittedParameters.value();
ACTS_VERBOSE("Fitted paramemeters for track " << itrack);
ACTS_VERBOSE(" position: " << params.position().transpose());
ACTS_VERBOSE(" momentum: " << params.momentum().transpose());
// Push the fitted parameters to the container
indexedParams.emplace(fitOutput.trackTip, std::move(params));
} else {
ACTS_DEBUG("No fitted paramemeters for track " << itrack);
}
// Create a SimMultiTrajectory
trajectories.emplace_back(std::move(fitOutput.fittedStates),
std::move(trackTips), std::move(indexedParams));
} else {
ACTS_WARNING("Fit failed for track " << itrack << " with error"
<< result.error());
// Fit failed, but still create an empty SimMultiTrajectory
trajectories.push_back(SimMultiTrajectory());
}
}
ctx.eventStore.add(m_cfg.outputTrajectories, std::move(trajectories));
return FW::ProcessCode::SUCCESS;
}
// This file is part of the Acts project.
//
// Copyright (C) 2019 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include <boost/program_options.hpp>
#include <iostream>
#include <map>
#include <random>
#include <stdexcept>
#include "ACTFW/Fitting/FittingAlgorithm.hpp"
#include "ACTFW/Plugins/BField/ScalableBField.hpp"
#include "Acts/Fitter/GainMatrixSmoother.hpp"
#include "Acts/Fitter/GainMatrixUpdater.hpp"
#include "Acts/Geometry/GeometryID.hpp"
#include "Acts/MagneticField/ConstantBField.hpp"
#include "Acts/MagneticField/InterpolatedBFieldMap.hpp"
#include "Acts/MagneticField/SharedBField.hpp"
#include "Acts/Propagator/EigenStepper.hpp"
#include "Acts/Propagator/Navigator.hpp"
#include "Acts/Propagator/Propagator.hpp"
#include "Acts/Surfaces/Surface.hpp"
#include "Acts/Utilities/Helpers.hpp"
#include "Acts/Utilities/ParameterDefinitions.hpp"
namespace {
template <typename Fitter>
struct FitterFunctionImpl {
Fitter fitter;
FitterFunctionImpl(Fitter&& f) : fitter(std::move(f)) {}
FW::FittingAlgorithm::FitterResult operator()(
const std::vector<FW::SimSourceLink>& sourceLinks,
const FW::TrackParameters& initialParameters,
const Acts::KalmanFitterOptions<Acts::VoidOutlierFinder>& options) const {
return fitter.fit(sourceLinks, initialParameters, options);
};
};
} // namespace
FW::FittingAlgorithm::FitterFunction FW::FittingAlgorithm::makeFitterFunction(
std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry,
Options::BFieldVariant magneticField, Acts::Logging::Level lvl) {
using Updater = Acts::GainMatrixUpdater<Acts::BoundParameters>;
using Smoother = Acts::GainMatrixSmoother<Acts::BoundParameters>;
// unpack the magnetic field variant and instantiate the corresponding fitter.
return std::visit(
[trackingGeometry, lvl](auto&& inputField) -> FitterFunction {
// each entry in the variant is already a shared_ptr
// need ::element_type to get the real magnetic field type
using InputMagneticField =
typename std::decay_t<decltype(inputField)>::element_type;
using MagneticField = Acts::SharedBField<InputMagneticField>;
using Stepper = Acts::EigenStepper<MagneticField>;
using Navigator = Acts::Navigator;
using Propagator = Acts::Propagator<Stepper, Navigator>;
using Fitter = Acts::KalmanFitter<Propagator, Updater, Smoother>;
// construct all components for the fitter
MagneticField field(std::move(inputField));
Stepper stepper(std::move(field));
Navigator navigator(trackingGeometry);
navigator.resolvePassive = false;
navigator.resolveMaterial = true;
navigator.resolveSensitive = true;
Propagator propagator(std::move(stepper), std::move(navigator));
Fitter fitter(std::move(propagator),
Acts::getDefaultLogger("KalmanFitter", lvl));
// build the fitter functions. owns the fitter object.
return FitterFunctionImpl<Fitter>(std::move(fitter));
},
std::move(magneticField));
}
add_library(
ActsExamplesGeant4 SHARED
src/EventAction.cpp
src/GdmlDetectorConstruction.cpp
src/GeantinoRecording.cpp
src/PrimaryGeneratorAction.cpp
src/RunAction.cpp
src/SteppingAction.cpp)
target_compile_definitions(
ActsExamplesGeant4
PUBLIC ${Geant4_DEFINITIONS})
target_include_directories(
ActsExamplesGeant4
SYSTEM PUBLIC ${Geant4_INCLUDE_DIRS})
target_include_directories(
ActsExamplesGeant4
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>)
target_link_libraries(
ActsExamplesGeant4
PUBLIC
ActsCore ActsExamplesDetectorDD4hep ActsExamplesFramework
${Geant4_LIBRARIES})
install(
TARGETS ActsExamplesGeant4
LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR})
// This file is part of the Acts project.
//
// Copyright (C) 2017-2020 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <G4VUserDetectorConstruction.hh>
#include <globals.hh>
namespace ActsExamples {
/// Construct the Geant4 detector from a Gdml file.
class GdmlDetectorConstruction final : public G4VUserDetectorConstruction {
public:
/// @param path is the path to the Gdml file
GdmlDetectorConstruction(std::string path);
/// Read the file and parse it to construct the Geant4 description.
G4VPhysicalVolume* Construct() final override;
private:
std::string m_path;
};
} // namespace ActsExamples
// This file is part of the Acts project.
//
// Copyright (C) 2017-2020 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <memory>
#include <mutex>
#include "ACTFW/Framework/BareAlgorithm.hpp"
#include "Acts/Propagator/MaterialInteractor.hpp"
#include "Acts/Utilities/Definitions.hpp"
#include "Acts/Utilities/Logger.hpp"
class G4RunManager;
class G4VUserDetectorConstruction;
namespace ActsExamples {
using RecordedMaterial = Acts::MaterialInteractor::result_type;
/// A material track with start position and momentum and recorded material.
using RecordedMaterialTrack =
std::pair<std::pair<Acts::Vector3D, Acts::Vector3D>, RecordedMaterial>;
/// Records the simulation geometry using geantinos.
///
/// This initiates the Geant4 simulation, and creates and writes out
/// the MaterialTrack entities which are needed for material mapping.
class GeantinoRecording final : public FW::BareAlgorithm {
public:
struct Config {
/// Output collection for the generated material tracks.
std::string outputMaterialTracks = "geant-material-tracks";
/// Detector construction object.
std::unique_ptr<G4VUserDetectorConstruction> detectorConstruction;
/// The number of tracks per event.
size_t tracksPerEvent = 0;
/// random number seed 1.
int seed1 = 12345;
/// random number seed 2.
int seed2 = 45678;
};
GeantinoRecording(Config&& cfg, Acts::Logging::Level lvl);
~GeantinoRecording();
FW::ProcessCode execute(const FW::AlgorithmContext& ctx) const final override;
private:
Config m_cfg;
std::unique_ptr<G4RunManager> m_runManager;
// has to be mutable; algorithm interface enforces object constness
mutable std::mutex m_runManagerLock;
};
} // namespace ActsExamples
// This file is part of the Acts project.
//
// Copyright (C) 2017-2020 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "EventAction.hpp"
#include <G4Event.hh>
#include <G4RunManager.hh>
#include <stdexcept>
#include "PrimaryGeneratorAction.hpp"
#include "SteppingAction.hpp"
using namespace ActsExamples;
EventAction* EventAction::s_instance = nullptr;
EventAction* EventAction::instance() {
return s_instance;
}
EventAction::EventAction() : G4UserEventAction() {
if (s_instance) {
throw std::logic_error("Attempted to duplicate the EventAction singleton");
} else {
s_instance = this;
}
}
EventAction::~EventAction() {
s_instance = nullptr;
}
void EventAction::BeginOfEventAction(const G4Event*) {
// reset the collection of material steps
SteppingAction::instance()->clear();
}
void EventAction::EndOfEventAction(const G4Event* event) {
const auto* rawPos = event->GetPrimaryVertex();
// access the initial direction of the track
G4ThreeVector rawDir = PrimaryGeneratorAction::instance()->direction();
// create the RecordedMaterialTrack
Acts::RecordedMaterialTrack mtrecord;
mtrecord.first.first =
Acts::Vector3D(rawPos->GetX0(), rawPos->GetY0(), rawPos->GetZ0());
mtrecord.first.second = Acts::Vector3D(rawDir.x(), rawDir.y(), rawDir.z());
mtrecord.second.materialInteractions =
SteppingAction::instance()->materialSteps();
// write out the RecordedMaterialTrack of one event
m_materialTracks.push_back(mtrecord);
}
/// Clear the recorded data.
void EventAction::clear() {
//
}
/// Access the recorded material tracks.
///
/// This only contains valid data after the end-of-event action has been
/// executed.
const std::vector<Acts::RecordedMaterialTrack>& EventAction::materialTracks()
const {
return m_materialTracks;
}
// This file is part of the Acts project.
//
// Copyright (C) 2017-2020 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#pragma once
#include <G4UserEventAction.hh>
#include <globals.hh>
#include <memory>
#include "ACTFW/EventData/SimHit.hpp"
#include "Acts/Propagator/MaterialInteractor.hpp"
namespace ActsExamples {
class SteppingAction;
/// @class EventAction
///
/// @brief Writes out material track records
///
/// The EventAction class is the realization of the Geant4 class
/// G4UserEventAction and is writing out the collected RecordedMaterialTrack
/// entities needed for material mapping once per event.
///
class EventAction final : public G4UserEventAction {
public:
/// Static access method
static EventAction* instance();
/// Construct the action and ensure singleton usage.
EventAction();
~EventAction() final override;
/// Interface method for begin of the event
/// @param event is the G4Event to be processed
/// @note resets the material step action
void BeginOfEventAction(const G4Event* event) final override;
/// Interface method for end of event
/// @param event is the G4Event to be processed
/// @note this method is writing out the material track records
void EndOfEventAction(const G4Event* event) final override;
/// Clear the recorded data.
void clear();
/// Access the recorded material tracks.
///
/// This only contains valid data after the end-of-event action has been
/// executed.
const std::vector<Acts::RecordedMaterialTrack>& materialTracks() const;
private:
/// Instance of the EventAction
static EventAction* s_instance;
/// The materialTrackWriter
std::vector<Acts::RecordedMaterialTrack> m_materialTracks;
};
} // namespace ActsExamples
// This file is part of the Acts project.
//
// Copyright (C) 2017-2020 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.
#include "ActsExamples/Geant4/GdmlDetectorConstruction.hpp"
#include <G4GDMLParser.hh>
using namespace ActsExamples;
GdmlDetectorConstruction::GdmlDetectorConstruction(std::string path)
: G4VUserDetectorConstruction(), m_path(std::move(path)) {}
G4VPhysicalVolume* GdmlDetectorConstruction::Construct() {
G4GDMLParser parser;
// TODO how to handle errors
parser.Read(m_path);
return parser.GetWorldVolume();
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment