From 59c1560740867d589a3c94f68ef3bb93d325f198 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 10:42:55 +0100 Subject: [PATCH 01/12] resovle merge conflicts --- Examples/Scripts/Python/sim_digi_odd.py | 234 ++++++++++++++++++++++++ 1 file changed, 234 insertions(+) create mode 100644 Examples/Scripts/Python/sim_digi_odd.py diff --git a/Examples/Scripts/Python/sim_digi_odd.py b/Examples/Scripts/Python/sim_digi_odd.py new file mode 100644 index 00000000000..3d06a93fa63 --- /dev/null +++ b/Examples/Scripts/Python/sim_digi_odd.py @@ -0,0 +1,234 @@ +#!/usr/bin/env python3 + +import os +import argparse +import pathlib + +import acts +import acts.examples +from acts.examples.simulation import ( + addParticleGun, + MomentumConfig, + EtaConfig, + PhiConfig, + ParticleConfig, + addPythia8, + addFatras, + addGeant4, + ParticleSelectorConfig, + addDigitization, + addParticleSelection, +) + +from acts.examples.odd import getOpenDataDetector, getOpenDataDetectorDirectory + +u = acts.UnitConstants + +parser = argparse.ArgumentParser(description="Sim-digi chain with the OpenDataDetector") +parser.add_argument( + "--output", + "-o", + help="Output directory", + type=pathlib.Path, + default=pathlib.Path.cwd() / "odd_output", +) +parser.add_argument("--events", "-n", help="Number of events", type=int, default=100) +parser.add_argument("--skip", "-s", help="Number of events", type=int, default=0) +parser.add_argument("--edm4hep", help="Use edm4hep inputs", type=pathlib.Path) +parser.add_argument( + "--geant4", help="Use Geant4 instead of fatras", action="store_true" +) +parser.add_argument( + "--ttbar", + help="Use Pythia8 (ttbar, pile-up 200) instead of particle gun", + action="store_true", +) +parser.add_argument( + "--ttbar-pu", + help="Number of pile-up events for ttbar", + type=int, + default=200, +) +parser.add_argument( + "--gun-multiplicity", + help="Multiplicity of the particle gun", + type=int, + default=200, +) +parser.add_argument( + "--gun-eta-range", + nargs="+", + help="Eta range of the particle gun", + type=float, + default=[-3.0, 3.0], +) +parser.add_argument( + "--gun-pt-range", + nargs="+", + help="Pt range of the particle gun (GeV)", + type=float, + default=[0.1 * u.GeV, 100 * u.GeV], +) +parser.add_argument( + "--rnd-seed", + help="Random seed", + type=int, + default=42, +) +parser.add_argument( + "--digi-config", + help="Digitization configuration file", + type=str, + default="", +) +args = parser.parse_args() + +outputDir = args.output +geoDir = getOpenDataDetectorDirectory() + +oddDigiConfig = args.digi_config + +detector, trackingGeometry, decorators = getOpenDataDetector( + odd_dir=geoDir, mdecorator=None +) +field = acts.ConstantBField(acts.Vector3(0.0, 0.0, 2.0 * u.T)) +rnd = acts.examples.RandomNumbers(seed=args.rnd_seed) + +s = acts.examples.Sequencer( + events=args.events, + skip=args.skip, + numThreads=1 if args.geant4 else -1, + outputDir=str(outputDir), +) + +if args.edm4hep: + import acts.examples.edm4hep + + edm4hepReader = acts.examples.edm4hep.EDM4hepReader( + inputPath=str(args.edm4hep), + inputSimHits=[ + "PixelBarrelReadout", + "PixelEndcapReadout", + "ShortStripBarrelReadout", + "ShortStripEndcapReadout", + "LongStripBarrelReadout", + "LongStripEndcapReadout", + ], + outputParticlesGenerator="particles_input", + outputParticlesInitial="particles_initial", + outputParticlesFinal="particles_final", + outputSimHits="simhits", + graphvizOutput="graphviz", + dd4hepDetector=detector, + trackingGeometry=trackingGeometry, + sortSimHitsInTime=True, + level=acts.logging.INFO, + ) + s.addReader(edm4hepReader) + s.addWhiteboardAlias("particles", edm4hepReader.config.outputParticlesGenerator) + + addParticleSelection( + s, + config=ParticleSelectorConfig( + rho=(0.0, 24 * u.mm), + absZ=(0.0, 1.0 * u.m), + eta=(-3.0, 3.0), + pt=(150 * u.MeV, None), + removeNeutral=True, + ), + inputParticles="particles", + outputParticles="particles_selected", + ) +else: + if not args.ttbar: + addParticleGun( + s, + MomentumConfig( + args.gun_pt_range[0] * u.GeV, + args.gun_pt_range[1] * u.GeV, + transverse=True, + ), + EtaConfig(args.gun_eta_range[0], args.gun_eta_range[1]), + PhiConfig(0.0, 360.0 * u.degree), + ParticleConfig(4, acts.PdgParticle.eMuon, randomizeCharge=True), + vtxGen=acts.examples.GaussianVertexGenerator( + mean=acts.Vector4(0, 0, 0, 0), + stddev=acts.Vector4( + 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 1.0 * u.ns + ), + ), + multiplicity=args.gun_multiplicity, + rnd=rnd, + ) + else: + addPythia8( + s, + hardProcess=["Top:qqbar2ttbar=on"], + npileup=args.ttbar_pu, + vtxGen=acts.examples.GaussianVertexGenerator( + mean=acts.Vector4(0, 0, 0, 0), + stddev=acts.Vector4( + 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 5.0 * u.ns + ), + ), + rnd=rnd, + outputDirRoot=outputDir, + # outputDirCsv=outputDir, + ) + + if args.geant4: + if s.config.numThreads != 1: + raise ValueError("Geant 4 simulation does not support multi-threading") + + # Pythia can sometime simulate particles outside the world volume, a cut on the Z of the track help mitigate this effect + # Older version of G4 might not work, this as has been tested on version `geant4-11-00-patch-03` + # For more detail see issue #1578 + addGeant4( + s, + detector, + trackingGeometry, + field, + preSelectParticles=ParticleSelectorConfig( + rho=(0.0, 24 * u.mm), + absZ=(0.0, 1.0 * u.m), + eta=(-3.0, 3.0), + pt=(150 * u.MeV, None), + removeNeutral=True, + ), + outputDirRoot=outputDir, + outputDirCsv=outputDir, + rnd=rnd, + killVolume=trackingGeometry.worldVolume, + killAfterTime=25 * u.ns, + ) + else: + addFatras( + s, + trackingGeometry, + field, + preSelectParticles=ParticleSelectorConfig( + rho=(0.0, 24 * u.mm), + absZ=(0.0, 1.0 * u.m), + eta=(-3.0, 3.0), + pt=(150 * u.MeV, None), + removeNeutral=True, + ) + if args.ttbar + else ParticleSelectorConfig(), + enableInteractions=True, + outputDirRoot=outputDir if args.output_root else None, + outputDirCsv=outputDir if args.output_csv else None, + rnd=rnd, + ) + +addDigitization( + s, + trackingGeometry, + field, + digiConfigFile=oddDigiConfig, + outputDirRoot=outputDir, + outputDirCsv=outputDir, + rnd=rnd, +) + +s.run() From 175def2d19d764b4179f189d8f93b51bd888237b Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 8 May 2024 22:02:30 +0200 Subject: [PATCH 02/12] add a few more options From b29d145cb59a8f400338768a24e4f3ff5354f83d Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 10:47:46 +0100 Subject: [PATCH 03/12] merge it --- Examples/Io/Obj/CMakeLists.txt | 1 + .../Obj/ObjPropagationStepsWriter.hpp | 0 .../ActsExamples/Io/Obj/ObjSimHitWriter.hpp | 77 +++++++++++++ .../Obj/ObjTrackingGeometryWriter.hpp | 0 .../{Plugins => Io}/Obj/ObjWriterOptions.hpp | 0 Examples/Io/Obj/src/ObjSimHitWriter.cpp | 106 ++++++++++++++++++ .../Io/Obj/src/ObjTrackingGeometryWriter.cpp | 2 +- .../Python/python/acts/examples/simulation.py | 24 +++- Examples/Python/src/Output.cpp | 5 + Examples/Scripts/Python/full_chain_odd.py | 8 ++ 10 files changed, 221 insertions(+), 2 deletions(-) rename Examples/Io/Obj/include/ActsExamples/{Plugins => Io}/Obj/ObjPropagationStepsWriter.hpp (100%) create mode 100644 Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp rename Examples/Io/Obj/include/ActsExamples/{Plugins => Io}/Obj/ObjTrackingGeometryWriter.hpp (100%) rename Examples/Io/Obj/include/ActsExamples/{Plugins => Io}/Obj/ObjWriterOptions.hpp (100%) create mode 100644 Examples/Io/Obj/src/ObjSimHitWriter.cpp diff --git a/Examples/Io/Obj/CMakeLists.txt b/Examples/Io/Obj/CMakeLists.txt index 9ca888f726a..bce7bf19e6b 100644 --- a/Examples/Io/Obj/CMakeLists.txt +++ b/Examples/Io/Obj/CMakeLists.txt @@ -2,6 +2,7 @@ add_library( ActsExamplesIoObj SHARED src/ObjTrackingGeometryWriter.cpp + src/ObjSimHitWriter.cpp src/ObjPropagationStepsWriter.cpp ) diff --git a/Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjPropagationStepsWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp similarity index 100% rename from Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjPropagationStepsWriter.hpp rename to Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp diff --git a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp new file mode 100644 index 00000000000..7a1da17a6f1 --- /dev/null +++ b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp @@ -0,0 +1,77 @@ +// This file is part of the Acts project. +// +// Copyright (C) 2024 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 "Acts/Utilities/Logger.hpp" +#include "ActsExamples/EventData/SimHit.hpp" +#include "ActsExamples/Framework/ProcessCode.hpp" +#include "ActsExamples/Framework/WriterT.hpp" + +#include +#include +#include + +namespace ActsExamples { +struct AlgorithmContext; + +/// Write out a simhit collection before detector digitization +/// as wavefront obj files. +/// +/// This writes one file per event containing information about the +/// global space points, momenta (before and after interaction) and hit index +/// into the configured output directory. By default it writes to the +/// current working directory. Files are named using the following schema +/// +/// event000000001-.obj +/// event000000002-.obj +/// ...class ObjSimHitWriter final : public WriterT { +class ObjSimHitWriter : public WriterT { + public: + struct Config { + /// Input sim hit collection to write. + std::string inputSimHits; + /// Where to place output files + std::string outputDir; + /// Output filename stem. + std::string outputStem = "simhits"; + /// Number of decimal digits for floating point precision in output. + std::size_t outputPrecision = std::numeric_limits::max_digits10; + /// Draw line connections between hits + bool drawConnections = true; + }; + + /// Construct the particle writer. + /// + /// @param config is the configuration object + /// @param level is the logging level + ObjSimHitWriter(const Config& config, Acts::Logging::Level level); + + /// Ensure underlying file is closed. + ~ObjSimHitWriter() override = default; + + /// End-of-run hook + ProcessCode finalize() override; + + /// Get readonly access to the config parameters + const Config& config() const { return m_cfg; } + + protected: + /// Type-specific write implementation. + /// + /// @param[in] ctx is the algorithm context + /// @param[in] hits are the hits to be written + ProcessCode writeT(const AlgorithmContext& ctx, + const SimHitContainer& hits) override; + + private: + Config m_cfg; + std::mutex m_writeMutex; +}; + +} // namespace ActsExamples diff --git a/Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjTrackingGeometryWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp similarity index 100% rename from Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjTrackingGeometryWriter.hpp rename to Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp diff --git a/Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjWriterOptions.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjWriterOptions.hpp similarity index 100% rename from Examples/Io/Obj/include/ActsExamples/Plugins/Obj/ObjWriterOptions.hpp rename to Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjWriterOptions.hpp diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp new file mode 100644 index 00000000000..d1d72c55c7f --- /dev/null +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -0,0 +1,106 @@ +// This file is part of the Acts project. +// +// Copyright (C) 2024 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/Io/Obj/ObjSimHitWriter.hpp" + +#include "Acts/Definitions/Algebra.hpp" +#include "Acts/Definitions/Common.hpp" +#include "Acts/Definitions/Units.hpp" +#include "Acts/Geometry/GeometryIdentifier.hpp" +#include "ActsExamples/EventData/SimHit.hpp" +#include "ActsExamples/Framework/AlgorithmContext.hpp" +#include "ActsExamples/Utilities/Paths.hpp" +#include "ActsFatras/EventData/Barcode.hpp" +#include "ActsFatras/EventData/Hit.hpp" + +#include +#include +#include +#include + +ActsExamples::ObjSimHitWriter::ObjSimHitWriter( + const ActsExamples::ObjSimHitWriter::Config& config, + Acts::Logging::Level level) + : WriterT(config.inputSimHits, "ObjSimHitWriter", level), m_cfg(config) { + // inputSimHits is already checked by base constructor + if (m_cfg.outputStem.empty()) { + throw std::invalid_argument("Missing output filename stem"); + } +} + +ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( + const AlgorithmContext& ctx, const ActsExamples::SimHitContainer& simHits) { + // ensure exclusive access to tree/file while writing + std::lock_guard lock(m_writeMutex); + + // open per-event file for all simhit components + std::string pathSimHit = perEventFilepath( + m_cfg.outputDir, m_cfg.outputStem + ".obj", ctx.eventNumber); + + std::ofstream os(pathSimHit, std::ofstream::out | std::ofstream::trunc); + if (!os) { + throw std::ios_base::failure("Could not open '" + pathSimHit + + "' to write"); + } + // Initialize the vertex counter + unsigned int vCounter = 0; + + if (!m_cfg.drawConnections) { + // Write data from internal immplementation + for (const auto& simHit : simHits) { + // local simhit information in global coord. + const Acts::Vector4& globalPos4 = simHit.fourPosition(); + + os << "v " << globalPos4[Acts::ePos0] / Acts::UnitConstants::mm << " " + << globalPos4[Acts::ePos1] / Acts::UnitConstants::mm << " " + << globalPos4[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; + + } // end simHit loop + } else { + // We need to associate first + std::unordered_map> particleHits; + // pre-loop over hits + for (const auto& simHit : simHits) { + if (particleHits.find(simHit.particleId().value()) == + particleHits.end()) { + particleHits[simHit.particleId().value()] = {}; + } + particleHits[simHit.particleId().value()].push_back( + simHit.fourPosition()); + } + // Draw loop + unsigned int lOffset = 1; + for (auto [pId, pHits] : particleHits) { + // Draw the particle hits + std::sort(pHits.begin(), pHits.end(), + [](const Acts::Vector4& a, const Acts::Vector4& b) { + return a[Acts::eTime] < b[Acts::eTime]; + }); + os << "o particle_" << pId << std::endl; + for (const auto& hit : pHits) { + os << "v " << hit[Acts::ePos0] / Acts::UnitConstants::mm << " " + << hit[Acts::ePos1] / Acts::UnitConstants::mm << " " + << hit[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; + } + // Draw the line + for (unsigned int iv = lOffset + 1; iv < lOffset + pHits.size(); ++iv) { + os << "l " << iv - 1 << " " << iv << '\n'; + } + os << '\n'; + // Increase the offset count + lOffset += pHits.size(); + } + } + os.close(); + + return ActsExamples::ProcessCode::SUCCESS; +} + +ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::finalize() { + return ActsExamples::ProcessCode::SUCCESS; +} diff --git a/Examples/Io/Obj/src/ObjTrackingGeometryWriter.cpp b/Examples/Io/Obj/src/ObjTrackingGeometryWriter.cpp index 4937757c85f..d2dbe52e745 100644 --- a/Examples/Io/Obj/src/ObjTrackingGeometryWriter.cpp +++ b/Examples/Io/Obj/src/ObjTrackingGeometryWriter.cpp @@ -6,7 +6,7 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. -#include "ActsExamples/Plugins/Obj/ObjTrackingGeometryWriter.hpp" +#include "ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp" #include "Acts/Utilities/Logger.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" diff --git a/Examples/Python/python/acts/examples/simulation.py b/Examples/Python/python/acts/examples/simulation.py index 65b6820d8b8..73a9d9b1f64 100644 --- a/Examples/Python/python/acts/examples/simulation.py +++ b/Examples/Python/python/acts/examples/simulation.py @@ -420,6 +420,7 @@ def addFatras( outputSimHits: str = "simhits", outputDirCsv: Optional[Union[Path, str]] = None, outputDirRoot: Optional[Union[Path, str]] = None, + outputDirObj: Optional[Union[Path, str]] = None, logLevel: Optional[acts.logging.Level] = None, ) -> None: """This function steers the detector simulation using Fatras @@ -444,6 +445,8 @@ def addFatras( the output folder for the Csv output, None triggers no output outputDirRoot : Path|str, path, None the output folder for the Root output, None triggers no output + outputDirObj : Path|str, path, None + the output folder for the Obj output, None triggers no output """ customLogLevel = acts.examples.defaultLogging(s, logLevel) @@ -507,6 +510,7 @@ def addFatras( particlesPostSelected, outputDirCsv, outputDirRoot, + outputDirObj, logLevel, ) @@ -519,6 +523,7 @@ def addSimWriters( particlesSimulated: str = "particles_simulated", outputDirCsv: Optional[Union[Path, str]] = None, outputDirRoot: Optional[Union[Path, str]] = None, + outputDirObj: Optional[Union[Path, str]] = None, logLevel: Optional[acts.logging.Level] = None, ) -> None: customLogLevel = acts.examples.defaultLogging(s, logLevel) @@ -563,6 +568,19 @@ def addSimWriters( ) ) + if outputDirObj is not None: + outputDirObj = Path(outputDirObj) + if not outputDirObj.exists(): + outputDirObj.mkdir() + s.addWriter( + acts.examples.ObjSimHitWriter( + level=customLogLevel(), + inputSimHits=simHits, + outputDir=str(outputDirObj), + outputStem="hits", + ) + ) + def getG4DetectorConstructionFactory( detector: Any, @@ -620,6 +638,7 @@ def addGeant4( keepParticlesWithoutHits=True, outputDirCsv: Optional[Union[Path, str]] = None, outputDirRoot: Optional[Union[Path, str]] = None, + outputDirObj: Optional[Union[Path, str]] = None, logLevel: Optional[acts.logging.Level] = None, killVolume: Optional[acts.Volume] = None, killAfterTime: float = float("inf"), @@ -647,6 +666,8 @@ def addGeant4( the output folder for the Csv output, None triggers no output outputDirRoot : Path|str, path, None the output folder for the Root output, None triggers no output + outputDirObj : Path|str, path, None + the output folder for the Obj output, None triggers no output killVolume: acts.Volume, None if given, particles are killed when going outside this volume. killAfterTime: float @@ -736,7 +757,8 @@ def addGeant4( particlesPostSelected, outputDirCsv, outputDirRoot, - logLevel, + outputDirObj, + logLevel=logLevel, ) return s diff --git a/Examples/Python/src/Output.cpp b/Examples/Python/src/Output.cpp index 4281f2db26e..6c92f05a6d0 100644 --- a/Examples/Python/src/Output.cpp +++ b/Examples/Python/src/Output.cpp @@ -47,6 +47,7 @@ #include "ActsExamples/Io/Root/VertexNTupleWriter.hpp" #include "ActsExamples/MaterialMapping/IMaterialWriter.hpp" #include "ActsExamples/Plugins/Obj/ObjPropagationStepsWriter.hpp" +#include "ActsExamples/Plugins/Obj/ObjSimHitWriter.hpp" #include "ActsExamples/Plugins/Obj/ObjTrackingGeometryWriter.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp" @@ -111,6 +112,10 @@ void addOutput(Context& ctx) { "ObjPropagationStepsWriter", collection, outputDir, outputScalor, outputPrecision); + ACTS_PYTHON_DECLARE_WRITER(ActsExamples::ObjSimHitWriter, mex, + "ObjSimHitWriter", inputSimHits, outputDir, + outputStem, outputPrecision, drawConnections); + { auto c = py::class_(m, "ViewConfig").def(py::init<>()); diff --git a/Examples/Scripts/Python/full_chain_odd.py b/Examples/Scripts/Python/full_chain_odd.py index 8b86402b820..7ddaa02054f 100755 --- a/Examples/Scripts/Python/full_chain_odd.py +++ b/Examples/Scripts/Python/full_chain_odd.py @@ -136,6 +136,12 @@ default=True, action=argparse.BooleanOptionalAction, ) +parser.add_argument( + "--output-obj", + help="Switch obj output on/off", + default=True, + action=argparse.BooleanOptionalAction, +) args = parser.parse_args() @@ -277,6 +283,7 @@ ), outputDirRoot=outputDir if args.output_root else None, outputDirCsv=outputDir if args.output_csv else None, + outputDirObj=outputDir if args.output_obj else None, rnd=rnd, killVolume=trackingGeometry.highestTrackingVolume, killAfterTime=25 * u.ns, @@ -305,6 +312,7 @@ enableInteractions=True, outputDirRoot=outputDir if args.output_root else None, outputDirCsv=outputDir if args.output_csv else None, + outputDirObj=outputDir if args.output_obj else None, rnd=rnd, ) From 66cea67614263a9075e57221ceb39110e6fa843c Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 10:48:41 +0100 Subject: [PATCH 04/12] merge conflict resolving --- Examples/Python/src/Output.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Examples/Python/src/Output.cpp b/Examples/Python/src/Output.cpp index 6c92f05a6d0..57725116db4 100644 --- a/Examples/Python/src/Output.cpp +++ b/Examples/Python/src/Output.cpp @@ -45,10 +45,10 @@ #include "ActsExamples/Io/Root/TrackFinderPerformanceWriter.hpp" #include "ActsExamples/Io/Root/TrackFitterPerformanceWriter.hpp" #include "ActsExamples/Io/Root/VertexNTupleWriter.hpp" +#include "ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp" +#include "ActsExamples/Io/Obj/ObjSimHitWriter.hpp" +#include "ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp" #include "ActsExamples/MaterialMapping/IMaterialWriter.hpp" -#include "ActsExamples/Plugins/Obj/ObjPropagationStepsWriter.hpp" -#include "ActsExamples/Plugins/Obj/ObjSimHitWriter.hpp" -#include "ActsExamples/Plugins/Obj/ObjTrackingGeometryWriter.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp" From e9409a79597c57be693270280e0cb71836efceba Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Thu, 9 May 2024 00:07:39 +0200 Subject: [PATCH 05/12] removing unused variable --- Examples/Io/Obj/src/ObjSimHitWriter.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index d1d72c55c7f..752e944d5ec 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -47,9 +47,8 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( throw std::ios_base::failure("Could not open '" + pathSimHit + "' to write"); } - // Initialize the vertex counter - unsigned int vCounter = 0; + // Only hit plotting if (!m_cfg.drawConnections) { // Write data from internal immplementation for (const auto& simHit : simHits) { From 5b5d3a5125fa17d59850175c884beef11808c179 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Mon, 13 May 2024 14:55:15 +0200 Subject: [PATCH 06/12] add spline fit to smoothen the curves --- .../ActsExamples/Io/Obj/ObjSimHitWriter.hpp | 9 ++ Examples/Io/Obj/src/ObjSimHitWriter.cpp | 109 +++++++++++++++--- 2 files changed, 99 insertions(+), 19 deletions(-) diff --git a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp index 7a1da17a6f1..d52a6787079 100644 --- a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp +++ b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp @@ -8,6 +8,7 @@ #pragma once +#include "Acts/Definitions/Units.hpp" #include "Acts/Utilities/Logger.hpp" #include "ActsExamples/EventData/SimHit.hpp" #include "ActsExamples/Framework/ProcessCode.hpp" @@ -17,6 +18,8 @@ #include #include +using namespace Acts::UnitLiterals; + namespace ActsExamples { struct AlgorithmContext; @@ -44,6 +47,12 @@ class ObjSimHitWriter : public WriterT { std::size_t outputPrecision = std::numeric_limits::max_digits10; /// Draw line connections between hits bool drawConnections = true; + /// Momentum threshold for hits + Acts::ActsScalar momentumThreshold = 0.05_GeV; + /// Momentum threshold for trajectories + Acts::ActsScalar momentumThresholdTraj = 0.05_GeV; + /// Number of points to interpolate between hits + unsigned int nInterpolatedPoints = 10; }; /// Construct the particle writer. diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index 752e944d5ec..126efee5ca7 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -10,7 +10,6 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/Common.hpp" -#include "Acts/Definitions/Units.hpp" #include "Acts/Geometry/GeometryIdentifier.hpp" #include "ActsExamples/EventData/SimHit.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" @@ -23,6 +22,33 @@ #include #include +#include + +namespace { + +template +std::vector interpolatedPoints( + const std::vector& inputs, unsigned int nPoints) { + Eigen::MatrixXd points(3, inputs.size()); + for (unsigned int i = 0; i < inputs.size(); ++i) { + points.col(i) = inputs[i].template head<3>().transpose(); + } + Eigen::Spline spline3D = + Eigen::SplineFitting>::Interpolate( + points, 2); + + std::vector output; + Acts::ActsScalar step = 1. / (nPoints - 1); + for (unsigned int i = 0; i < nPoints; ++i) { + Acts::ActsScalar t = i * step; + Eigen::Vector3d point = spline3D(t); + output.push_back(Acts::Vector3(point[0], point[1], point[2])); + } + return output; +} + +} // namespace + ActsExamples::ObjSimHitWriter::ObjSimHitWriter( const ActsExamples::ObjSimHitWriter::Config& config, Acts::Logging::Level level) @@ -41,11 +67,16 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( // open per-event file for all simhit components std::string pathSimHit = perEventFilepath( m_cfg.outputDir, m_cfg.outputStem + ".obj", ctx.eventNumber); + std::string pathSimTrajectory = perEventFilepath( + m_cfg.outputDir, m_cfg.outputStem + "_trajectory.obj", ctx.eventNumber); - std::ofstream os(pathSimHit, std::ofstream::out | std::ofstream::trunc); - if (!os) { - throw std::ios_base::failure("Could not open '" + pathSimHit + - "' to write"); + std::ofstream osHits(pathSimHit, std::ofstream::out | std::ofstream::trunc); + std::ofstream osTrajectory(pathSimTrajectory, + std::ofstream::out | std::ofstream::trunc); + + if (!osHits || !osTrajectory) { + throw std::ios_base::failure("Could not open '" + pathSimHit + "' or '" + + pathSimTrajectory + "' to write"); } // Only hit plotting @@ -55,16 +86,34 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( // local simhit information in global coord. const Acts::Vector4& globalPos4 = simHit.fourPosition(); - os << "v " << globalPos4[Acts::ePos0] / Acts::UnitConstants::mm << " " - << globalPos4[Acts::ePos1] / Acts::UnitConstants::mm << " " - << globalPos4[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; + osHits << "v " << globalPos4[Acts::ePos0] / Acts::UnitConstants::mm << " " + << globalPos4[Acts::ePos1] / Acts::UnitConstants::mm << " " + << globalPos4[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; } // end simHit loop } else { // We need to associate first std::unordered_map> particleHits; - // pre-loop over hits + // Pre-loop over hits ... write those below threshold for (const auto& simHit : simHits) { + Acts::ActsScalar momentum = simHit.momentum4Before().head<3>().norm(); + if (momentum < m_cfg.momentumThreshold) { + ACTS_VERBOSE("Skipping : Hit below threshold: " << momentum); + continue; + } else if (momentum < m_cfg.momentumThresholdTraj) { + ACTS_VERBOSE( + "Skipping (trajectory): Hit below threshold: " << momentum); + osHits << "v " + << simHit.fourPosition()[Acts::ePos0] / Acts::UnitConstants::mm + << " " + << simHit.fourPosition()[Acts::ePos1] / Acts::UnitConstants::mm + << " " + << simHit.fourPosition()[Acts::ePos2] / Acts::UnitConstants::mm + << std::endl; + continue; + } + ACTS_VERBOSE("Accepting : Hit above threshold: " << momentum); + if (particleHits.find(simHit.particleId().value()) == particleHits.end()) { particleHits[simHit.particleId().value()] = {}; @@ -74,28 +123,50 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( } // Draw loop unsigned int lOffset = 1; - for (auto [pId, pHits] : particleHits) { + for (auto& [pId, pHits] : particleHits) { // Draw the particle hits std::sort(pHits.begin(), pHits.end(), [](const Acts::Vector4& a, const Acts::Vector4& b) { return a[Acts::eTime] < b[Acts::eTime]; }); - os << "o particle_" << pId << std::endl; + + osHits << "o particle_" << pId << std::endl; for (const auto& hit : pHits) { - os << "v " << hit[Acts::ePos0] / Acts::UnitConstants::mm << " " - << hit[Acts::ePos1] / Acts::UnitConstants::mm << " " - << hit[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; + osHits << "v " << hit[Acts::ePos0] / Acts::UnitConstants::mm << " " + << hit[Acts::ePos1] / Acts::UnitConstants::mm << " " + << hit[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; + } + osHits << '\n'; + + // Interpolate the points + std::vector trajectory; + if (pHits.size() < 3) { + for (const auto& hit : pHits) { + trajectory.push_back(hit.template head<3>()); + } + } else { + trajectory = + interpolatedPoints(pHits, pHits.size() * m_cfg.nInterpolatedPoints); + } + + osTrajectory << "o particle_trajectory_" << pId << std::endl; + for (const auto& hit : trajectory) { + osTrajectory << "v " << hit[Acts::ePos0] / Acts::UnitConstants::mm + << " " << hit[Acts::ePos1] / Acts::UnitConstants::mm << " " + << hit[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; } // Draw the line - for (unsigned int iv = lOffset + 1; iv < lOffset + pHits.size(); ++iv) { - os << "l " << iv - 1 << " " << iv << '\n'; + for (unsigned int iv = lOffset + 1; iv < lOffset + trajectory.size(); + ++iv) { + osTrajectory << "l " << iv - 1 << " " << iv << '\n'; } - os << '\n'; + osTrajectory << '\n'; // Increase the offset count - lOffset += pHits.size(); + lOffset += trajectory.size(); } } - os.close(); + osHits.close(); + osTrajectory.close(); return ActsExamples::ProcessCode::SUCCESS; } From 44e59d34eebc6a125fbfc8b333ae6176871b0775 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 11:05:04 +0100 Subject: [PATCH 07/12] rebase from May --- .../ActsExamples/Io/Obj/ObjSimHitWriter.hpp | 12 +++++++++-- .../Io/Obj/src/ObjPropagationStepsWriter.cpp | 2 +- Examples/Io/Obj/src/ObjSimHitWriter.cpp | 21 ++++++++++++------- Examples/Python/src/Output.cpp | 6 +++--- Examples/Scripts/Python/sim_digi_odd.py | 20 ++++++++++-------- 5 files changed, 39 insertions(+), 22 deletions(-) diff --git a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp index d52a6787079..852991800ad 100644 --- a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp +++ b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp @@ -1,3 +1,11 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/. + // This file is part of the Acts project. // // Copyright (C) 2024 CERN for the benefit of the Acts project @@ -48,9 +56,9 @@ class ObjSimHitWriter : public WriterT { /// Draw line connections between hits bool drawConnections = true; /// Momentum threshold for hits - Acts::ActsScalar momentumThreshold = 0.05_GeV; + double momentumThreshold = 0.05_GeV; /// Momentum threshold for trajectories - Acts::ActsScalar momentumThresholdTraj = 0.05_GeV; + double momentumThresholdTraj = 0.05_GeV; /// Number of points to interpolate between hits unsigned int nInterpolatedPoints = 10; }; diff --git a/Examples/Io/Obj/src/ObjPropagationStepsWriter.cpp b/Examples/Io/Obj/src/ObjPropagationStepsWriter.cpp index 6d0388984cc..0f196ce162c 100644 --- a/Examples/Io/Obj/src/ObjPropagationStepsWriter.cpp +++ b/Examples/Io/Obj/src/ObjPropagationStepsWriter.cpp @@ -6,7 +6,7 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. -#include "ActsExamples/Plugins/Obj/ObjPropagationStepsWriter.hpp" +#include "ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp" namespace ActsExamples { diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index 126efee5ca7..5a3a38954c2 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -1,3 +1,11 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/. + // This file is part of the Acts project. // // Copyright (C) 2024 CERN for the benefit of the Acts project @@ -33,14 +41,13 @@ std::vector interpolatedPoints( for (unsigned int i = 0; i < inputs.size(); ++i) { points.col(i) = inputs[i].template head<3>().transpose(); } - Eigen::Spline spline3D = - Eigen::SplineFitting>::Interpolate( - points, 2); + Eigen::Spline spline3D = + Eigen::SplineFitting>::Interpolate(points, 2); std::vector output; - Acts::ActsScalar step = 1. / (nPoints - 1); + double step = 1. / (nPoints - 1); for (unsigned int i = 0; i < nPoints; ++i) { - Acts::ActsScalar t = i * step; + double t = i * step; Eigen::Vector3d point = spline3D(t); output.push_back(Acts::Vector3(point[0], point[1], point[2])); } @@ -93,10 +100,10 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( } // end simHit loop } else { // We need to associate first - std::unordered_map> particleHits; + std::unordered_map> particleHits; // Pre-loop over hits ... write those below threshold for (const auto& simHit : simHits) { - Acts::ActsScalar momentum = simHit.momentum4Before().head<3>().norm(); + double momentum = simHit.momentum4Before().head<3>().norm(); if (momentum < m_cfg.momentumThreshold) { ACTS_VERBOSE("Skipping : Hit below threshold: " << momentum); continue; diff --git a/Examples/Python/src/Output.cpp b/Examples/Python/src/Output.cpp index 57725116db4..251ca2b6357 100644 --- a/Examples/Python/src/Output.cpp +++ b/Examples/Python/src/Output.cpp @@ -25,6 +25,9 @@ #include "ActsExamples/Io/Csv/CsvTrackParameterWriter.hpp" #include "ActsExamples/Io/Csv/CsvTrackWriter.hpp" #include "ActsExamples/Io/Csv/CsvTrackingGeometryWriter.hpp" +#include "ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp" +#include "ActsExamples/Io/Obj/ObjSimHitWriter.hpp" +#include "ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp" #include "ActsExamples/Io/Root/RootBFieldWriter.hpp" #include "ActsExamples/Io/Root/RootMaterialTrackWriter.hpp" #include "ActsExamples/Io/Root/RootMaterialWriter.hpp" @@ -45,9 +48,6 @@ #include "ActsExamples/Io/Root/TrackFinderPerformanceWriter.hpp" #include "ActsExamples/Io/Root/TrackFitterPerformanceWriter.hpp" #include "ActsExamples/Io/Root/VertexNTupleWriter.hpp" -#include "ActsExamples/Io/Obj/ObjPropagationStepsWriter.hpp" -#include "ActsExamples/Io/Obj/ObjSimHitWriter.hpp" -#include "ActsExamples/Io/Obj/ObjTrackingGeometryWriter.hpp" #include "ActsExamples/MaterialMapping/IMaterialWriter.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupReader.hpp" #include "ActsExamples/TrackFinding/ITrackParamsLookupWriter.hpp" diff --git a/Examples/Scripts/Python/sim_digi_odd.py b/Examples/Scripts/Python/sim_digi_odd.py index 3d06a93fa63..24f03cd9da3 100644 --- a/Examples/Scripts/Python/sim_digi_odd.py +++ b/Examples/Scripts/Python/sim_digi_odd.py @@ -206,15 +206,17 @@ s, trackingGeometry, field, - preSelectParticles=ParticleSelectorConfig( - rho=(0.0, 24 * u.mm), - absZ=(0.0, 1.0 * u.m), - eta=(-3.0, 3.0), - pt=(150 * u.MeV, None), - removeNeutral=True, - ) - if args.ttbar - else ParticleSelectorConfig(), + preSelectParticles=( + ParticleSelectorConfig( + rho=(0.0, 24 * u.mm), + absZ=(0.0, 1.0 * u.m), + eta=(-3.0, 3.0), + pt=(150 * u.MeV, None), + removeNeutral=True, + ) + if args.ttbar + else ParticleSelectorConfig() + ), enableInteractions=True, outputDirRoot=outputDir if args.output_root else None, outputDirCsv=outputDir if args.output_csv else None, From ba77cbb6fc5232a72421838461c440f2123abb4b Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 11:13:43 +0100 Subject: [PATCH 08/12] remove sim_digi_odd.py --- Examples/Scripts/Python/sim_digi_odd.py | 236 ------------------------ 1 file changed, 236 deletions(-) delete mode 100644 Examples/Scripts/Python/sim_digi_odd.py diff --git a/Examples/Scripts/Python/sim_digi_odd.py b/Examples/Scripts/Python/sim_digi_odd.py deleted file mode 100644 index 24f03cd9da3..00000000000 --- a/Examples/Scripts/Python/sim_digi_odd.py +++ /dev/null @@ -1,236 +0,0 @@ -#!/usr/bin/env python3 - -import os -import argparse -import pathlib - -import acts -import acts.examples -from acts.examples.simulation import ( - addParticleGun, - MomentumConfig, - EtaConfig, - PhiConfig, - ParticleConfig, - addPythia8, - addFatras, - addGeant4, - ParticleSelectorConfig, - addDigitization, - addParticleSelection, -) - -from acts.examples.odd import getOpenDataDetector, getOpenDataDetectorDirectory - -u = acts.UnitConstants - -parser = argparse.ArgumentParser(description="Sim-digi chain with the OpenDataDetector") -parser.add_argument( - "--output", - "-o", - help="Output directory", - type=pathlib.Path, - default=pathlib.Path.cwd() / "odd_output", -) -parser.add_argument("--events", "-n", help="Number of events", type=int, default=100) -parser.add_argument("--skip", "-s", help="Number of events", type=int, default=0) -parser.add_argument("--edm4hep", help="Use edm4hep inputs", type=pathlib.Path) -parser.add_argument( - "--geant4", help="Use Geant4 instead of fatras", action="store_true" -) -parser.add_argument( - "--ttbar", - help="Use Pythia8 (ttbar, pile-up 200) instead of particle gun", - action="store_true", -) -parser.add_argument( - "--ttbar-pu", - help="Number of pile-up events for ttbar", - type=int, - default=200, -) -parser.add_argument( - "--gun-multiplicity", - help="Multiplicity of the particle gun", - type=int, - default=200, -) -parser.add_argument( - "--gun-eta-range", - nargs="+", - help="Eta range of the particle gun", - type=float, - default=[-3.0, 3.0], -) -parser.add_argument( - "--gun-pt-range", - nargs="+", - help="Pt range of the particle gun (GeV)", - type=float, - default=[0.1 * u.GeV, 100 * u.GeV], -) -parser.add_argument( - "--rnd-seed", - help="Random seed", - type=int, - default=42, -) -parser.add_argument( - "--digi-config", - help="Digitization configuration file", - type=str, - default="", -) -args = parser.parse_args() - -outputDir = args.output -geoDir = getOpenDataDetectorDirectory() - -oddDigiConfig = args.digi_config - -detector, trackingGeometry, decorators = getOpenDataDetector( - odd_dir=geoDir, mdecorator=None -) -field = acts.ConstantBField(acts.Vector3(0.0, 0.0, 2.0 * u.T)) -rnd = acts.examples.RandomNumbers(seed=args.rnd_seed) - -s = acts.examples.Sequencer( - events=args.events, - skip=args.skip, - numThreads=1 if args.geant4 else -1, - outputDir=str(outputDir), -) - -if args.edm4hep: - import acts.examples.edm4hep - - edm4hepReader = acts.examples.edm4hep.EDM4hepReader( - inputPath=str(args.edm4hep), - inputSimHits=[ - "PixelBarrelReadout", - "PixelEndcapReadout", - "ShortStripBarrelReadout", - "ShortStripEndcapReadout", - "LongStripBarrelReadout", - "LongStripEndcapReadout", - ], - outputParticlesGenerator="particles_input", - outputParticlesInitial="particles_initial", - outputParticlesFinal="particles_final", - outputSimHits="simhits", - graphvizOutput="graphviz", - dd4hepDetector=detector, - trackingGeometry=trackingGeometry, - sortSimHitsInTime=True, - level=acts.logging.INFO, - ) - s.addReader(edm4hepReader) - s.addWhiteboardAlias("particles", edm4hepReader.config.outputParticlesGenerator) - - addParticleSelection( - s, - config=ParticleSelectorConfig( - rho=(0.0, 24 * u.mm), - absZ=(0.0, 1.0 * u.m), - eta=(-3.0, 3.0), - pt=(150 * u.MeV, None), - removeNeutral=True, - ), - inputParticles="particles", - outputParticles="particles_selected", - ) -else: - if not args.ttbar: - addParticleGun( - s, - MomentumConfig( - args.gun_pt_range[0] * u.GeV, - args.gun_pt_range[1] * u.GeV, - transverse=True, - ), - EtaConfig(args.gun_eta_range[0], args.gun_eta_range[1]), - PhiConfig(0.0, 360.0 * u.degree), - ParticleConfig(4, acts.PdgParticle.eMuon, randomizeCharge=True), - vtxGen=acts.examples.GaussianVertexGenerator( - mean=acts.Vector4(0, 0, 0, 0), - stddev=acts.Vector4( - 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 1.0 * u.ns - ), - ), - multiplicity=args.gun_multiplicity, - rnd=rnd, - ) - else: - addPythia8( - s, - hardProcess=["Top:qqbar2ttbar=on"], - npileup=args.ttbar_pu, - vtxGen=acts.examples.GaussianVertexGenerator( - mean=acts.Vector4(0, 0, 0, 0), - stddev=acts.Vector4( - 0.0125 * u.mm, 0.0125 * u.mm, 55.5 * u.mm, 5.0 * u.ns - ), - ), - rnd=rnd, - outputDirRoot=outputDir, - # outputDirCsv=outputDir, - ) - - if args.geant4: - if s.config.numThreads != 1: - raise ValueError("Geant 4 simulation does not support multi-threading") - - # Pythia can sometime simulate particles outside the world volume, a cut on the Z of the track help mitigate this effect - # Older version of G4 might not work, this as has been tested on version `geant4-11-00-patch-03` - # For more detail see issue #1578 - addGeant4( - s, - detector, - trackingGeometry, - field, - preSelectParticles=ParticleSelectorConfig( - rho=(0.0, 24 * u.mm), - absZ=(0.0, 1.0 * u.m), - eta=(-3.0, 3.0), - pt=(150 * u.MeV, None), - removeNeutral=True, - ), - outputDirRoot=outputDir, - outputDirCsv=outputDir, - rnd=rnd, - killVolume=trackingGeometry.worldVolume, - killAfterTime=25 * u.ns, - ) - else: - addFatras( - s, - trackingGeometry, - field, - preSelectParticles=( - ParticleSelectorConfig( - rho=(0.0, 24 * u.mm), - absZ=(0.0, 1.0 * u.m), - eta=(-3.0, 3.0), - pt=(150 * u.MeV, None), - removeNeutral=True, - ) - if args.ttbar - else ParticleSelectorConfig() - ), - enableInteractions=True, - outputDirRoot=outputDir if args.output_root else None, - outputDirCsv=outputDir if args.output_csv else None, - rnd=rnd, - ) - -addDigitization( - s, - trackingGeometry, - field, - digiConfigFile=oddDigiConfig, - outputDirRoot=outputDir, - outputDirCsv=outputDir, - rnd=rnd, -) - -s.run() From f549616a8fdc9bca3b43df27c90862c143700831 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 11:18:24 +0100 Subject: [PATCH 09/12] adding output obj --- Examples/Scripts/Python/geant4.py | 1 + 1 file changed, 1 insertion(+) diff --git a/Examples/Scripts/Python/geant4.py b/Examples/Scripts/Python/geant4.py index 42cb5db95bb..65e154c2704 100755 --- a/Examples/Scripts/Python/geant4.py +++ b/Examples/Scripts/Python/geant4.py @@ -36,6 +36,7 @@ def runGeant4( field, outputDirCsv=outputDir / "csv", outputDirRoot=outputDir, + outputDirObj=outputDir / "obj", rnd=rnd, materialMappings=materialMappings, volumeMappings=volumeMappings, From e903ba2c888b1b4639ad84219ed163348c0d1507 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 13:10:30 +0100 Subject: [PATCH 10/12] make SonarCloud happy --- .../Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp | 6 ++---- Examples/Io/Obj/src/ObjSimHitWriter.cpp | 2 +- 2 files changed, 3 insertions(+), 5 deletions(-) diff --git a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp index 852991800ad..fc9b3772e06 100644 --- a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp +++ b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp @@ -26,8 +26,6 @@ #include #include -using namespace Acts::UnitLiterals; - namespace ActsExamples { struct AlgorithmContext; @@ -56,9 +54,9 @@ class ObjSimHitWriter : public WriterT { /// Draw line connections between hits bool drawConnections = true; /// Momentum threshold for hits - double momentumThreshold = 0.05_GeV; + double momentumThreshold = 0.05 * Acts::UnitConstants::GeV; /// Momentum threshold for trajectories - double momentumThresholdTraj = 0.05_GeV; + double momentumThresholdTraj = 0.05 * Acts::UnitConstants::GeV; /// Number of points to interpolate between hits unsigned int nInterpolatedPoints = 10; }; diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index 5a3a38954c2..6c7b4ff107d 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -69,7 +69,7 @@ ActsExamples::ObjSimHitWriter::ObjSimHitWriter( ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( const AlgorithmContext& ctx, const ActsExamples::SimHitContainer& simHits) { // ensure exclusive access to tree/file while writing - std::lock_guard lock(m_writeMutex); + std::scoped_lock lock(m_writeMutex); // open per-event file for all simhit components std::string pathSimHit = perEventFilepath( From c184d19b32058f2ec29a08077b96aab667b27646 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 16:10:20 +0100 Subject: [PATCH 11/12] addressing PR comments --- .../ActsExamples/Io/Obj/ObjSimHitWriter.hpp | 24 +++----- Examples/Io/Obj/src/ObjSimHitWriter.cpp | 57 +++++++++++-------- 2 files changed, 40 insertions(+), 41 deletions(-) diff --git a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp index fc9b3772e06..8ccd92645b9 100644 --- a/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp +++ b/Examples/Io/Obj/include/ActsExamples/Io/Obj/ObjSimHitWriter.hpp @@ -6,14 +6,6 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. -// This file is part of the Acts project. -// -// Copyright (C) 2024 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 "Acts/Definitions/Units.hpp" @@ -29,17 +21,17 @@ namespace ActsExamples { struct AlgorithmContext; -/// Write out a simhit collection before detector digitization -/// as wavefront obj files. +/// Write out a simhit collection before detector digitization as wavefront obj +/// file(s per event). /// -/// This writes one file per event containing information about the -/// global space points, momenta (before and after interaction) and hit index -/// into the configured output directory. By default it writes to the -/// current working directory. Files are named using the following schema +/// This writes two files per event, one for the hits one for the trajectory. +/// The latter can be smoothed using a spline interpolation. /// /// event000000001-.obj +/// event000000001-_trajectory.obj /// event000000002-.obj -/// ...class ObjSimHitWriter final : public WriterT { +/// event000000002-_trajectory.obj +/// class ObjSimHitWriter : public WriterT { public: struct Config { @@ -58,7 +50,7 @@ class ObjSimHitWriter : public WriterT { /// Momentum threshold for trajectories double momentumThresholdTraj = 0.05 * Acts::UnitConstants::GeV; /// Number of points to interpolate between hits - unsigned int nInterpolatedPoints = 10; + std::size_t nInterpolatedPoints = 10; }; /// Construct the particle writer. diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index 6c7b4ff107d..377d0795e2c 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -6,14 +6,6 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. -// This file is part of the Acts project. -// -// Copyright (C) 2024 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/Io/Obj/ObjSimHitWriter.hpp" #include "Acts/Definitions/Algebra.hpp" @@ -34,22 +26,37 @@ namespace { +/// @brief Helper function to interpolate points +/// +/// @tparam input_vector_type +/// @param inputs input vector points +/// @param nPoints number of interpolation points +/// +/// @return std::vector interpolated points template std::vector interpolatedPoints( - const std::vector& inputs, unsigned int nPoints) { - Eigen::MatrixXd points(3, inputs.size()); - for (unsigned int i = 0; i < inputs.size(); ++i) { - points.col(i) = inputs[i].template head<3>().transpose(); - } - Eigen::Spline spline3D = - Eigen::SplineFitting>::Interpolate(points, 2); - + const std::vector& inputs, std::size_t nPoints) { std::vector output; - double step = 1. / (nPoints - 1); - for (unsigned int i = 0; i < nPoints; ++i) { - double t = i * step; - Eigen::Vector3d point = spline3D(t); - output.push_back(Acts::Vector3(point[0], point[1], point[2])); + + if (nPoints < 1) { + // No interpolation done return simply the output vector + for (const auto& input : inputs) { + output.push_back(input.template head<3>()); + } + + } else { + Eigen::MatrixXd points(3, inputs.size()); + for (std::size_t i = 0; i < inputs.size(); ++i) { + points.col(i) = inputs[i].template head<3>().transpose(); + } + Eigen::Spline spline3D = + Eigen::SplineFitting>::Interpolate(points, 2); + + double step = 1. / (nPoints - 1); + for (std::size_t i = 0; i < nPoints; ++i) { + double t = i * step; + output.push_back(spline3D(t)); + } } return output; } @@ -105,7 +112,7 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( for (const auto& simHit : simHits) { double momentum = simHit.momentum4Before().head<3>().norm(); if (momentum < m_cfg.momentumThreshold) { - ACTS_VERBOSE("Skipping : Hit below threshold: " << momentum); + ACTS_VERBOSE("Skipping: Hit below threshold: " << momentum); continue; } else if (momentum < m_cfg.momentumThresholdTraj) { ACTS_VERBOSE( @@ -119,7 +126,7 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( << std::endl; continue; } - ACTS_VERBOSE("Accepting : Hit above threshold: " << momentum); + ACTS_VERBOSE("Accepting: Hit above threshold: " << momentum); if (particleHits.find(simHit.particleId().value()) == particleHits.end()) { @@ -129,7 +136,7 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( simHit.fourPosition()); } // Draw loop - unsigned int lOffset = 1; + std::size_t lOffset = 1; for (auto& [pId, pHits] : particleHits) { // Draw the particle hits std::sort(pHits.begin(), pHits.end(), @@ -163,7 +170,7 @@ ActsExamples::ProcessCode ActsExamples::ObjSimHitWriter::writeT( << hit[Acts::ePos2] / Acts::UnitConstants::mm << std::endl; } // Draw the line - for (unsigned int iv = lOffset + 1; iv < lOffset + trajectory.size(); + for (std::size_t iv = lOffset + 1; iv < lOffset + trajectory.size(); ++iv) { osTrajectory << "l " << iv - 1 << " " << iv << '\n'; } From 8ded30c7ff25a2aa50d4ec4807c07ee0e5c9dc80 Mon Sep 17 00:00:00 2001 From: Andreas Salzburger Date: Wed, 4 Dec 2024 16:34:24 +0100 Subject: [PATCH 12/12] Update Examples/Io/Obj/src/ObjSimHitWriter.cpp Co-authored-by: Alexander J. Pfleger <70842573+AJPfleger@users.noreply.github.com> --- Examples/Io/Obj/src/ObjSimHitWriter.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Examples/Io/Obj/src/ObjSimHitWriter.cpp b/Examples/Io/Obj/src/ObjSimHitWriter.cpp index 377d0795e2c..564ff9c5409 100644 --- a/Examples/Io/Obj/src/ObjSimHitWriter.cpp +++ b/Examples/Io/Obj/src/ObjSimHitWriter.cpp @@ -38,7 +38,7 @@ std::vector interpolatedPoints( const std::vector& inputs, std::size_t nPoints) { std::vector output; - if (nPoints < 1) { + if (nPoints < 2) { // No interpolation done return simply the output vector for (const auto& input : inputs) { output.push_back(input.template head<3>());