Skip to content

Commit 19d966b

Browse files
committed
Make per-track random seeding work for Geant4
SimCutParams.trackSeed was so far a no-op with Geant4 because seeding only occurred in the Geant3 stack-pop path. Add the missing seeding hook in O2MCApplicationBase::PreTrack(), which Geant4 invokes for both primary and secondary tracks on their first step. Geant3 continues to seed at stack-pop time. This preserves existing behaviour and avoids moving seeding to a later hook, which measurably reduces reproducibility. The PreTrack seed is derived from the engine track state rather than Stack::GetCurrentTrack(), which is only reliable for Geant4 primaries. Also add diagnostics to detect missing seed propagation by counting successful seed applications and warning at the end of the event if track seeding was requested but never performed. Verified on both Geant3 and Geant4: per-track seeding now works as intended and existing Geant3 behaviour remains unchanged.
1 parent 70803ce commit 19d966b

4 files changed

Lines changed: 94 additions & 2 deletions

File tree

Detectors/Base/include/DetectorsBase/VMCSeederService.h

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -35,13 +35,17 @@ class VMCSeederService
3535

3636
void setSeed() const; // will propagate seed to the VMC engines
3737

38+
/// how often a seed was propagated; lets callers detect a silent no-op
39+
unsigned long long getSeedCount() const { return mSeedCount; }
40+
3841
typedef std::function<void()> SeederFcn;
3942

4043
private:
4144
VMCSeederService();
4245
void initSeederFunction(TVirtualMC const*);
4346

44-
SeederFcn mSeederFcn; // the just-in-time compiled function talking to the VMC engines
47+
SeederFcn mSeederFcn; // the just-in-time compiled function talking to the VMC engines
48+
mutable unsigned long long mSeedCount{0}; // number of setSeed() calls
4549
};
4650

4751
} // namespace base

Detectors/Base/src/VMCSeederService.cxx

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,4 +50,5 @@ void VMCSeederService::setSeed() const
5050
// This is ok since in any case gRandom->SetSeed(seed); gRandom->GetSeed() != seed;
5151
gRandom->Rndm();
5252
mSeederFcn();
53+
++mSeedCount;
5354
}

Steer/include/Steer/O2MCApplicationBase.h

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -68,6 +68,11 @@ class O2MCApplicationBase : public FairMCApplication
6868
// keeping track of volumeIds and volume names
6969

7070
double mLongestTrackTime = 0;
71+
bool mTrackSeedWarned{false}; // whether we already complained that seeding never fired
72+
73+
/// whether this engine needs per-track seeding in PreTrack (Geant3 seeds at
74+
/// stack-pop time instead, see O2MCApplicationBase::seedsInPreTrack)
75+
bool seedsInPreTrack() const;
7176
/// some common parts of finishEvent
7277
void finishEventCommon();
7378
TrackRefFcn mTrackRefFcn; // a function hook that gets (optionally) called during Stepping

Steer/src/O2MCApplication.cxx

Lines changed: 83 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,10 @@
4646
#include <DetectorsBase/O2Tessellated.h>
4747
#include <unordered_set>
4848
#include "SimConfig/G4Params.h"
49+
#include "DetectorsBase/VMCSeederService.h" // per-track seeding of the engine
50+
#include <TLorentzVector.h>
51+
#include <TRandom.h>
52+
#include <cstring>
4953

5054
namespace o2
5155
{
@@ -119,9 +123,76 @@ void O2MCApplicationBase::Stepping()
119123
FairMCApplication::Stepping();
120124
}
121125

126+
namespace
127+
{
128+
// Hash of a track's initial state (vertex, global time, momentum, PDG). Used as
129+
// the random seed for that track, so that a track's random stream depends only
130+
// on the track itself and not on how many randoms earlier tracks happened to
131+
// consume.
132+
//
133+
// The values are read from the transport engine, not from
134+
// o2::data::Stack::GetCurrentTrack(): under Geant4 the stack's "current track"
135+
// is only meaningful for primaries -- Stack::SetCurrentTrack() falls back to
136+
// mCurrentParticle0 (the last particle *pushed*) for anything beyond the
137+
// primary array, so every secondary would hash the wrong particle. Both engines
138+
// have the track's initial state loaded by the time PreTrack is called (Geant4
139+
// sets the step to kVertex first; Geant3 calls GLTRAC before GUTRAK).
140+
ULong_t hashCurrentTrack(TVirtualMC* vmc)
141+
{
142+
auto asLong = [](double x) {
143+
ULong_t l;
144+
std::memcpy(&l, &x, sizeof(l));
145+
return l;
146+
};
147+
148+
TLorentzVector pos, mom;
149+
vmc->TrackPosition(pos);
150+
vmc->TrackMomentum(mom);
151+
152+
ULong_t hash = asLong(pos.X());
153+
hash ^= asLong(pos.Y());
154+
hash ^= asLong(pos.Z());
155+
hash ^= asLong(pos.T());
156+
hash ^= asLong(mom.Px());
157+
hash ^= asLong(mom.Py());
158+
hash ^= asLong(mom.Pz());
159+
hash += (ULong_t)vmc->TrackPid();
160+
return hash;
161+
}
162+
} // namespace
163+
164+
bool O2MCApplicationBase::seedsInPreTrack() const
165+
{
166+
// Geant3 seeds at stack-pop time, in o2::data::Stack::PopNextTrack(). That is
167+
// strictly earlier than its PreTrack hook (gutrak) and measurably stronger:
168+
// with the TOF module removed from an otherwise identical setup, pop-time
169+
// seeding keeps all 603 ITS hits bit-identical, PreTrack seeding only 68 %.
170+
// Do not seed Geant3 here as well -- it is already covered, and reseeding a
171+
// second time mid-track would undo the first.
172+
static const bool inPreTrack = [this]() {
173+
const char* name = (fMC != nullptr) ? fMC->GetName() : "";
174+
return strncmp(name, "TGeant3", 7) != 0;
175+
}();
176+
return inPreTrack;
177+
}
178+
122179
void O2MCApplicationBase::PreTrack()
123180
{
124-
// dispatch first to function in FairRoot
181+
if (mCutParams.trackSeed && seedsInPreTrack()) {
182+
// Per-track seeding for engines that do not go through
183+
// o2::data::Stack::PopNextTrack(). Geant4 is one: it takes primaries via
184+
// PopPrimaryForTracking and keeps secondaries internally, so the stack hook
185+
// never fires and this is the only per-track hook available. It is called
186+
// for primaries and secondaries alike
187+
// (TG4TrackingAction::PreUserTrackingAction), and only on a track's first
188+
// step, so a suspended track is not reseeded mid-flight.
189+
auto hash = hashCurrentTrack(fMC);
190+
// TRandom::SetSeed(0) means "seed from the clock" -- never let that happen.
191+
gRandom->SetSeed(hash == 0 ? 1 : hash);
192+
o2::base::VMCSeederService::instance().setSeed();
193+
}
194+
195+
// dispatch now to function in FairRoot
125196
FairMCApplication::PreTrack();
126197
}
127198

@@ -309,6 +380,17 @@ void O2MCApplicationBase::finishEventCommon()
309380
header->setDetId2HitBitLUT(o2::base::Detector::getDetId2HitBitIndex());
310381

311382
static_cast<o2::data::Stack*>(GetStack())->updateEventStats();
383+
384+
// Per-track seeding used to be wired to a stack callback that one of the two
385+
// engines never invoked, and it failed silently. Never again: if it was asked
386+
// for and nothing was seeded, say so.
387+
if (mCutParams.trackSeed && o2::base::VMCSeederService::instance().getSeedCount() == 0 &&
388+
!mTrackSeedWarned) {
389+
mTrackSeedWarned = true;
390+
LOG(warn) << "Per-track seeding (SimCutParams.trackSeed) was requested but not a single track "
391+
"was seeded -- neither the stack nor the PreTrack hook fired for this engine. "
392+
"Seeding is NOT active.";
393+
}
312394
}
313395

314396
void O2MCApplicationBase::FinishEvent()

0 commit comments

Comments
 (0)