Skip to content

Commit

Permalink
feat: Add python vertex example + test (#1028)
Browse files Browse the repository at this point in the history
This test tries to model all the vertexing combinations, where I noticed that some of them don't seem to produce outputs that make sense to me. It might be that some combinations are not supported / never tested, and I added a warning to the example script if such a combination is configured.

Right now this is still WIP because it's rebased on #1027 since it uses the truth tracking for the test

@baschlag could you help and take a look at this to see if this implementation makes sense?
  • Loading branch information
paulgessinger committed Jan 10, 2022
1 parent 8b0005c commit 5e7485b
Show file tree
Hide file tree
Showing 7 changed files with 345 additions and 9 deletions.
1 change: 1 addition & 0 deletions Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.ipp
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ auto Acts::AdaptiveMultiVertexFinder<vfitter_t, sfinder_t>::find(
const VertexingOptions<InputTrack_t>& vertexingOptions,
State& /*state*/) const -> Result<std::vector<Vertex<InputTrack_t>>> {
if (allTracks.empty()) {
ACTS_ERROR("Empty track collection handed to find method");
return VertexingError::EmptyInput;
}
// Original tracks
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -117,16 +117,26 @@ ActsExamples::AdaptiveMultiVertexFinderAlgorithm::execute(
using VertexingOptions = Acts::VertexingOptions<Acts::BoundTrackParameters>;
VertexingOptions finderOpts(ctx.geoContext, ctx.magFieldContext);

// find vertices and measure elapsed time
std::vector<Acts::Vertex<Fitter::InputTrack_t>> vertices;

auto t1 = std::chrono::high_resolution_clock::now();
auto result = finder.find(inputTrackPointers, finderOpts, state);
auto t2 = std::chrono::high_resolution_clock::now();

if (not result.ok()) {
ACTS_ERROR("Error in vertex finder: " << result.error().message());
return ProcessCode::ABORT;
if (inputTrackParameters.empty()) {
ACTS_DEBUG("Empty track parameter collection found, skipping vertexing");
} else {
ACTS_DEBUG("Have " << inputTrackParameters.size()
<< " input track parameters, running vertexing");
// find vertices and measure elapsed time
auto result = finder.find(inputTrackPointers, finderOpts, state);

if (not result.ok()) {
ACTS_ERROR("Error in vertex finder: " << result.error().message());
return ProcessCode::ABORT;
}
vertices = *result;
}
auto vertices = *result;

auto t2 = std::chrono::high_resolution_clock::now();

// show some debug output
ACTS_INFO("Found " << vertices.size() << " vertices in event");
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ ActsExamples::ProcessCode ActsExamples::VertexFitterAlgorithm::execute(
for (const auto& protoVertex : protoVertices) {
// un-constrained fit requires at least two tracks
if ((not m_cfg.doConstrainedFit) and (protoVertex.size() < 2)) {
ACTS_WARNING(
ACTS_INFO(
"Skip un-constrained vertex fit on proto-vertex with less than two "
"tracks");
continue;
Expand Down
2 changes: 1 addition & 1 deletion Examples/Io/Root/src/RootVertexPerformanceWriter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -201,7 +201,7 @@ ActsExamples::ProcessCode ActsExamples::RootVertexPerformanceWriter::writeT(
"Total number of reconstructed tracks : " << inputFittedTracks.size());

if (associatedTruthParticles.size() != inputFittedTracks.size()) {
ACTS_WARNING(
ACTS_INFO(
"Number of fitted tracks and associated truth particles do not match. "
"Not able to match fitted tracks at reconstructed vertex to truth "
"vertex.");
Expand Down
6 changes: 6 additions & 0 deletions Examples/Python/tests/root_file_hashes.txt
Original file line number Diff line number Diff line change
Expand Up @@ -49,3 +49,9 @@ test_truth_tracking[1000.0]__performance_track_finder.root: 76a990d595b6e097da2b

test_digitization_example_input__measurements.root: ccc92f0ad538d1b62d98f19f947970bcc491843e54d8ffeed16ad2e226b8caee
test_digitization_example_input__particles.root: 78a89f365177423d0834ea6f1bd8afe1488e72b12a25066a20bd9050f5407860

test_vertex_fitting_reading[Truth-False-100]__performance_vertexing.root: 489b058ee78b4d01075fc2c594ffeed4af07c306122e577857b1b7e9357888cf
test_vertex_fitting_reading[Iterative-False-100]__performance_vertexing.root: 85754514dff6640401af5cafc8dffb7ffa02ffd9358f658d89c9c39b623111e3
test_vertex_fitting_reading[Iterative-True-100]__performance_vertexing.root: 5a9e2d202015fe0e16b266e55896162c94654f55f1ace95b295ba3aecfc6cab8
test_vertex_fitting_reading[AMVF-False-100]__performance_vertexing.root: 8fb55f5aceae330e7738986030d9d1726b9f0b825449d22f7a13f7580890b7d5
test_vertex_fitting_reading[AMVF-True-100]__performance_vertexing.root: 907fffdc3534a29c82afb83325bd6107c7b7f9a6114f4dd8bef5745e5729d993
111 changes: 111 additions & 0 deletions Examples/Python/tests/test_examples.py
Original file line number Diff line number Diff line change
Expand Up @@ -838,3 +838,114 @@ def test_ckf_tracks_example_truth_smeared(tmp_path, assert_root_hash):

assert len([f for f in csv.iterdir() if f.name.endswith("CKFtracks.csv")]) == events
assert all([f.stat().st_size > 300 for f in csv.iterdir()])


@pytest.mark.skipif(not dd4hepEnabled, reason="DD4hep not set up")
@pytest.mark.slow
# @pytest.mark.filterwarnings("ignore::UserWarning")
def test_vertex_fitting(tmp_path):
detector, trackingGeometry, decorators = getOpenDataDetector()

field = acts.ConstantBField(acts.Vector3(0, 0, 2 * u.T))

from vertex_fitting import runVertexFitting, VertexFinder

s = Sequencer(events=100)

runVertexFitting(
field,
vertexFinder=VertexFinder.Truth,
outputDir=Path.cwd(),
s=s,
)

alg = AssertCollectionExistsAlg(["fittedVertices"], name="check_alg")
s.addAlgorithm(alg)

s.run()
assert alg.events_seen == s.config.events


import itertools


@pytest.mark.parametrize(
"finder,inputTracks,entries",
[
("Truth", False, 100),
# ("Truth", True, 0), # this combination seems to be not working
("Iterative", False, 100),
("Iterative", True, 100),
("AMVF", False, 100),
("AMVF", True, 100),
],
)
@pytest.mark.filterwarnings("ignore::UserWarning")
def test_vertex_fitting_reading(
tmp_path, ptcl_gun, rng, finder, inputTracks, entries, assert_root_hash
):

ptcl_file = tmp_path / "particles.root"

detector, trackingGeometry, decorators = GenericDetector.create()
field = acts.ConstantBField(acts.Vector3(0, 0, 2 * u.T))

from vertex_fitting import runVertexFitting, VertexFinder

inputTrackSummary = None
if inputTracks:
from truth_tracking import runTruthTracking

s2 = Sequencer(numThreads=1, events=100)
runTruthTracking(
trackingGeometry,
field,
digiConfigFile=Path(
"Examples/Algorithms/Digitization/share/default-smearing-config-generic.json"
),
outputDir=tmp_path,
s=s2,
)
s2.run()
del s2
inputTrackSummary = tmp_path / "tracksummary_fitter.root"
assert inputTrackSummary.exists()
assert ptcl_file.exists()
else:
s0 = Sequencer(events=100, numThreads=1)
evGen = ptcl_gun(s0)
s0.addWriter(
RootParticleWriter(
level=acts.logging.INFO,
inputParticles=evGen.config.outputParticles,
filePath=str(ptcl_file),
)
)
s0.run()
del s0

assert ptcl_file.exists()

finder = VertexFinder[finder]

s3 = Sequencer(numThreads=1)

runVertexFitting(
field,
inputParticlePath=ptcl_file,
inputTrackSummary=inputTrackSummary,
outputDir=tmp_path,
vertexFinder=finder,
s=s3,
)

alg = AssertCollectionExistsAlg(["fittedVertices"], name="check_alg")
s3.addAlgorithm(alg)

s3.run()

vertexing_file = tmp_path / "performance_vertexing.root"
assert vertexing_file.exists()

assert_entries(vertexing_file, "vertexing", entries)
assert_root_hash(vertexing_file.name, vertexing_file)
208 changes: 208 additions & 0 deletions Examples/Scripts/Python/vertex_fitting.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,208 @@
#!/usr/bin/env python3
from pathlib import Path
from typing import Optional
import enum
import warnings

from acts.examples import (
Sequencer,
ParticleSelector,
ParticleSmearing,
TruthVertexFinder,
VertexFitterAlgorithm,
IterativeVertexFinderAlgorithm,
RootParticleReader,
AdaptiveMultiVertexFinderAlgorithm,
RootVertexPerformanceWriter,
RootTrajectorySummaryReader,
TrackSelector,
GenericDetector,
)

import acts

from acts import UnitConstants as u

from common import addPythia8


class VertexFinder(enum.Enum):
Truth = (1,)
AMVF = (2,)
Iterative = (3,)


def runVertexFitting(
field,
outputDir: Path,
outputRoot: bool = True,
inputParticlePath: Optional[Path] = None,
inputTrackSummary: Path = None,
vertexFinder: VertexFinder = VertexFinder.Truth,
s=None,
):
s = s or Sequencer(events=100, numThreads=-1)

logger = acts.logging.getLogger("VertexFittingExample")

rnd = acts.examples.RandomNumbers(seed=42)

if inputParticlePath is None:
logger.info("Generating particles using Pythia8")
evGen = addPythia8(s, rnd)
inputParticles = evGen.config.outputParticles
else:
logger.info("Reading particles from %s", inputParticlePath.resolve())
assert inputParticlePath.exists()
inputParticles = "particles_read"
s.addReader(
RootParticleReader(
level=acts.logging.INFO,
filePath=str(inputParticlePath.resolve()),
particleCollection=inputParticles,
orderedEvents=False,
)
)

selectedParticles = "particles_selected"
ptclSelector = ParticleSelector(
level=acts.logging.INFO,
inputParticles=inputParticles,
outputParticles=selectedParticles,
removeNeutral=True,
absEtaMax=2.5,
rhoMax=4.0 * u.mm,
ptMin=500 * u.MeV,
)
s.addAlgorithm(ptclSelector)

trackParameters = "trackparameters"
if inputTrackSummary is None or inputParticlePath is None:
logger.info("Using smeared particles")

ptclSmearing = ParticleSmearing(
level=acts.logging.INFO,
inputParticles=selectedParticles,
outputTrackParameters=trackParameters,
randomNumbers=rnd,
)
s.addAlgorithm(ptclSmearing)
associatedParticles = selectedParticles
else:
logger.info("Reading track summary from %s", inputTrackSummary.resolve())
assert inputTrackSummary.exists()
associatedParticles = "associatedTruthParticles"
trackSummaryReader = RootTrajectorySummaryReader(
level=acts.logging.VERBOSE,
outputTracks="fittedTrackParameters",
outputParticles=associatedParticles,
filePath=str(inputTrackSummary.resolve()),
orderedEvents=False,
)
s.addReader(trackSummaryReader)

s.addAlgorithm(
TrackSelector(
level=acts.logging.INFO,
inputTrackParameters=trackSummaryReader.config.outputTracks,
outputTrackParameters=trackParameters,
outputTrackIndices="outputTrackIndices",
removeNeutral=True,
absEtaMax=2.5,
loc0Max=4.0 * u.mm, # rho max
ptMin=500 * u.MeV,
)
)

logger.info("Using vertex finder: %s", vertexFinder.name)

outputVertices = "fittedVertices"
outputTime = ""
if vertexFinder == VertexFinder.Truth:
findVertices = TruthVertexFinder(
level=acts.logging.VERBOSE,
inputParticles=selectedParticles,
outputProtoVertices="protovertices",
excludeSecondaries=True,
)
s.addAlgorithm(findVertices)
fitVertices = VertexFitterAlgorithm(
level=acts.logging.VERBOSE,
bField=field,
inputTrackParameters=trackParameters,
inputProtoVertices=findVertices.config.outputProtoVertices,
outputVertices=outputVertices,
)
s.addAlgorithm(fitVertices)

elif vertexFinder == VertexFinder.Iterative:
findVertices = IterativeVertexFinderAlgorithm(
level=acts.logging.INFO,
bField=field,
inputTrackParameters=trackParameters,
outputProtoVertices="protovertices",
outputVertices=outputVertices,
)
s.addAlgorithm(findVertices)
elif vertexFinder == VertexFinder.AMVF:
outputTime = "outputTime"
findVertices = AdaptiveMultiVertexFinderAlgorithm(
level=acts.logging.INFO,
bField=field,
inputTrackParameters=trackParameters,
outputProtoVertices="protovertices",
outputVertices=outputVertices,
outputTime=outputTime,
)

s.addAlgorithm(findVertices)
else:
raise RuntimeError("Invalid finder argument")

if outputRoot:
if inputTrackSummary is None:
warnings.warn(
"Using inputTrackSummary == None with outputRoot: "
"This combination is not necessarily supported. "
"Please get in touch with us"
)
s.addWriter(
RootVertexPerformanceWriter(
level=acts.logging.INFO,
inputAllTruthParticles=inputParticles,
inputSelectedTruthParticles=selectedParticles,
inputAssociatedTruthParticles=associatedParticles,
inputFittedTracks=trackParameters,
inputVertices=outputVertices,
inputTime=outputTime,
treeName="vertexing",
filePath=str(outputDir / "performance_vertexing.root"),
)
)

return s


if "__main__" == __name__:
detector, trackingGeometry, decorators = GenericDetector.create()

field = acts.ConstantBField(acts.Vector3(0, 0, 2 * u.T))

inputParticlePath = Path("particles.root")
if not inputParticlePath.exists():
inputParticlePath = None

inputTrackSummary = None
for p in ("tracksummary_fitter.root", "tracksummary_ckf.root"):
p = Path(p)
if p.exists():
inputTrackSummary = p
break

runVertexFitting(
field,
vertexFinder=VertexFinder.Truth,
inputParticlePath=inputParticlePath,
inputTrackSummary=inputTrackSummary,
outputDir=Path.cwd(),
).run()

0 comments on commit 5e7485b

Please sign in to comment.