Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

feat: persist ACTS barcodes in G4 sim #1842

Merged
merged 10 commits into from
Feb 9, 2023
Original file line number Diff line number Diff line change
Expand Up @@ -42,6 +42,12 @@ class EventStoreRegistry {
std::unordered_map<size_t, Acts::RecordedMaterialTrack> materialTracks;
/// Particle hit count (for hit indexing)
std::unordered_map<SimBarcode, std::size_t> particleHitCount;
/// Track ID to Barcode mapping
std::unordered_map<unsigned int, SimBarcode> trackIdMapping;
/// Track ID to root Track ID mapping
std::unordered_map<unsigned int, unsigned int> trackIdRootId;
/// Track ID generation counter
std::unordered_map<unsigned int, unsigned int> trackIdGenerationCount;
};

EventStoreRegistry() = delete;
Expand Down
23 changes: 10 additions & 13 deletions Examples/Algorithms/Geant4/src/Geant4Simulation.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -168,33 +168,30 @@ ActsExamples::ProcessCode ActsExamples::Geant4Simulation::execute(
if (not m_cfg.outputParticlesInitial.empty() and
not m_cfg.outputParticlesFinal.empty()) {
// Initial state of particles
SimParticleContainer outputParticlesInitial;
outputParticlesInitial.insert(eventData.particlesInitial.begin(),
eventData.particlesInitial.end());
// Register to the event store
ctx.eventStore.add(m_cfg.outputParticlesInitial,
std::move(outputParticlesInitial));
SimParticleContainer(eventData.particlesInitial.begin(),
eventData.particlesInitial.end()));
// Final state of particles
SimParticleContainer outputParticlesFinal;
outputParticlesFinal.insert(eventData.particlesFinal.begin(),
eventData.particlesFinal.end());
// Register to the event store
ctx.eventStore.add(m_cfg.outputParticlesFinal,
std::move(outputParticlesFinal));
SimParticleContainer(eventData.particlesFinal.begin(),
eventData.particlesFinal.end()));
}

// Output handling: Simulated hits
if (not m_cfg.outputSimHits.empty()) {
SimHitContainer simHits;
simHits.insert(eventData.hits.begin(), eventData.hits.end());
// Register to the event store
ctx.eventStore.add(m_cfg.outputSimHits, std::move(simHits));
ctx.eventStore.add(
m_cfg.outputSimHits,
SimHitContainer(eventData.hits.begin(), eventData.hits.end()));
}

// Output handling: Material tracks
if (not m_cfg.outputMaterialTracks.empty()) {
ctx.eventStore.add(m_cfg.outputMaterialTracks,
std::move(eventData.materialTracks));
ctx.eventStore.add(
m_cfg.outputMaterialTracks,
decltype(eventData.materialTracks)(eventData.materialTracks));
}

return ActsExamples::ProcessCode::SUCCESS;
Expand Down
26 changes: 25 additions & 1 deletion Examples/Algorithms/Geant4/src/ParticleTrackingAction.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,8 @@ void ActsExamples::ParticleTrackingAction::PostUserTrackingAction(

ActsExamples::SimParticle ActsExamples::ParticleTrackingAction::convert(
const G4Track& aTrack) const {
auto& eventData = EventStoreRegistry::eventData();

// Unit conversions G4->::ACTS
constexpr double convertTime = Acts::UnitConstants::s / CLHEP::s;
constexpr double convertLength = Acts::UnitConstants::mm / CLHEP::mm;
Expand All @@ -47,12 +49,34 @@ ActsExamples::SimParticle ActsExamples::ParticleTrackingAction::convert(
G4double charge = particleDef->GetPDGCharge();
G4int pdg = particleDef->GetPDGEncoding();
G4int id = aTrack.GetTrackID();
G4int parentId = aTrack.GetParentID();
G4ThreeVector pPosition = convertLength * aTrack.GetPosition();
G4double pTime = convertTime * aTrack.GetGlobalTime();
G4ThreeVector pDirection = aTrack.GetMomentumDirection();
G4double p = convertEnergy * aTrack.GetKineticEnergy();

if (parentId != 0) {
eventData.trackIdRootId[id] = eventData.trackIdRootId[parentId];
} else {
eventData.trackIdRootId[id] = id;
}

SimBarcode particleId;
if (eventData.trackIdMapping.find(id) != eventData.trackIdMapping.end()) {
particleId = eventData.trackIdMapping[id];
} else {
if (eventData.trackIdRootId.find(id) != eventData.trackIdRootId.end()) {
auto rootId = eventData.trackIdRootId[id];
particleId = eventData.trackIdMapping[rootId];
particleId.setGeneration(++eventData.trackIdGenerationCount[rootId]);
eventData.trackIdMapping[id] = particleId;
} else {
ACTS_WARNING("could not find parent " << parentId << " of " << id);
}
}

// Now create the Particle
ActsExamples::SimParticle aParticle(SimBarcode(id), Acts::PdgParticle(pdg),
ActsExamples::SimParticle aParticle(particleId, Acts::PdgParticle(pdg),
charge, mass);
aParticle.setPosition4(pPosition[0], pPosition[1], pPosition[2], pTime);
aParticle.setDirection(pDirection[0], pDirection[1], pDirection[2]);
Expand Down
9 changes: 5 additions & 4 deletions Examples/Algorithms/Geant4/src/SensitiveSteppingAction.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,9 @@ void ActsExamples::SensitiveSteppingAction::UserSteppingAction(
constexpr double convertLength = Acts::UnitConstants::mm / CLHEP::mm;
constexpr double convertEnergy = Acts::UnitConstants::GeV / CLHEP::GeV;

// Retrieve the event data registry
auto& eventData = EventStoreRegistry::eventData();

// The particle after the step
G4Track* track = step->GetTrack();
G4PrimaryParticle* primaryParticle =
Expand Down Expand Up @@ -101,14 +104,12 @@ void ActsExamples::SensitiveSteppingAction::UserSteppingAction(
Acts::GeometryIdentifier::Value sGeoVal = std::stoul(volumeName);
Acts::GeometryIdentifier geoID(sGeoVal);

SimBarcode particleID(track->GetTrackID());
auto particleID = eventData.trackIdMapping[track->GetTrackID()];

Acts::Vector4 particlePosition(hX, hY, hZ, hT);
Acts::Vector4 beforeMomentum(mXpre, mYpre, mZpre, mEpre);
Acts::Vector4 afterMomentum(mXpost, mYpost, mZpost, mEpost);

// Retrieve the event data registry
auto& eventData = EventStoreRegistry::eventData();

// Increase counter (starts at 1 because of ++)
++eventData.particleHitCount[particleID];

Expand Down
19 changes: 13 additions & 6 deletions Examples/Algorithms/Geant4/src/SimParticleTranslation.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -65,30 +65,33 @@ void ActsExamples::SimParticleTranslation::GeneratePrimaries(G4Event* anEvent) {
G4PrimaryVertex* pVertex = nullptr;

// We are looping through the particles and flush per vertex
std::unique_ptr<SimParticle::Vector4> lastVertex = nullptr;
std::optional<Acts::Vector4> lastVertex;

constexpr double convertLength = CLHEP::mm / Acts::UnitConstants::mm;
constexpr double convertEnergy = CLHEP::GeV / Acts::UnitConstants::GeV;

unsigned int pCounter = 0;
unsigned int trackId = 1;
// Loop over the input partilces and run
for (const auto& part : inputParticles) {
auto currentVertex = part.fourPosition();
if (lastVertex == nullptr or not currentVertex.isApprox(*lastVertex)) {
if (not lastVertex or not currentVertex.isApprox(*lastVertex)) {
// Add the vertex to the event
if (pVertex != nullptr) {
anEvent->AddPrimaryVertex(pVertex);
ACTS_DEBUG("Flushing " << pCounter
<< " particles associated with vertex "
<< Acts::toString(*lastVertex));
<< lastVertex->transpose());
pCounter = 0;
}
lastVertex = std::make_unique<SimParticle::Vector4>(currentVertex);
lastVertex = currentVertex;
pVertex = new G4PrimaryVertex(
currentVertex[0] * convertLength, currentVertex[1] * convertLength,
currentVertex[2] * convertLength, currentVertex[3]);
}

// Add a new primary to the vertex

Acts::Vector4 mom4 = part.fourMomentum() * convertEnergy;

// Particle properties, may be forced to specific value
Expand Down Expand Up @@ -132,15 +135,19 @@ void ActsExamples::SimParticleTranslation::GeneratePrimaries(G4Event* anEvent) {
particle->SetMass(particleMass);
particle->Set4Momentum(mom4[0], mom4[1], mom4[2], mom4[3]);
particle->SetCharge(part.charge());
particle->SetTrackID(pCounter);
particle->SetTrackID(trackId++);

// Add the primary to the vertex
pVertex->SetPrimary(particle);

eventData.trackIdMapping[particle->GetTrackID()] = part.particleId();

++pCounter;
}
// Final vertex to be added
if (pVertex != nullptr) {
anEvent->AddPrimaryVertex(pVertex);
ACTS_DEBUG("Flushing " << pCounter << " particles associated with vertex "
<< Acts::toString(*lastVertex));
<< lastVertex->transpose());
}
}
6 changes: 3 additions & 3 deletions Examples/Python/tests/root_file_hashes.txt
Original file line number Diff line number Diff line change
Expand Up @@ -2,9 +2,9 @@ test_pythia8__pythia8_particles.root: 7bae8c6353d792ea2d54401334e2e57e486ec154ee
test_fatras__fatras_particles_final.root: 7ce243ae0e64b13d1f14c7b40f7d6a3bf9edd799f15bf0745b0b08c16ff11bf4
test_fatras__fatras_particles_initial.root: c72be956b187ea67cb5efe8fed6f7cdd6657631d25aa68105bcfe5af5baecc41
test_fatras__hits.root: 2f25d996715c6ed05e65a3f36df559de52e0a6c8ab9641fedd810474eaf5db11
test_geant4__fatras_particles_final.root: 944eed8b749e48031cd8332f7a6994880bd2f381b93d10e887e1e86387c3e5c3
test_geant4__fatras_particles_initial.root: 161525bb3c62cf9746f5b0528e2d65e9f6d32a4e015a52e2e6d8c6feb13d7753
test_geant4__hits.root: 3949a5c201c874c1366a09473c024986dc231baf8279b202107785c4fffb0347
test_geant4__fatras_particles_final.root: b6254065a9bf03fe538d8b39615bfb3936a5e3eef28c4ad67d59514b170d4f87
test_geant4__fatras_particles_initial.root: 077546f5598fbf32ea15dc1b2a4741c7bfe51d6d3b96421ff04b999aaf21c2fd
test_geant4__hits.root: 507461900225433b2452c66db9a9712bfe431f157bff9f98115b30d2b0ac361e
test_seeding__estimatedparams.root: 3d181aeda75f4a030ab015c90291c21448c61e513ad1e7bfcaf56317b7af1a68
test_seeding__performance_seeding_trees.root: b53253e5549c09282d68cf11716ff53fa68ea340e25277571b4bc52d43804469
test_seeding__particles.root: 8c22c826bf7b46de59f479b00e1aa54871570ba85a09d0b221527499f84c4be3
Expand Down