Skip to content

Commit 06083b1

Browse files
sawenzelsawenzel
committed
MCTrack: Generalize hit bit property
So far, every possible sensitive detector for ALICE Run3 had a dedicated bit assigned inside MCTracks. The number of such bit was fixed to be 22 maximally. In this way we could store information if a particular Monte Carlo track caused energy deposits in detectors. Since more and more detectors are being added to the code (ALICE3, ITS3, R&D), we'll soon have more than 22 DetIDs in our code which will conflict with the space available in MCTrack. Since enlarging the MCTrack class in memory is to be avoided, this commit proposes a solution based on an intermediate lookup table which dynamically maps the detector ID to a particular bit. The lookup table is setup during geometry construction and stored in MCEventHeaders for offline use. This solution should scale much further than before. The limit is now 22 simultaneous sensitive detectors in the geometry and not 22 total sensitive detectors in the code.
1 parent 769b8fc commit 06083b1

12 files changed

Lines changed: 145 additions & 71 deletions

File tree

‎DataFormats/simulation/include/SimulationDataFormat/MCEventHeader.h‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -86,6 +86,9 @@ class MCEventHeader : public FairMCEventHeader
8686

8787
MCEventStats& getMCEventStats() { return mEventStats; }
8888

89+
void setDetId2HitBitLUT(std::vector<int> const& v) { mDetId2HitBitIndex = v; }
90+
std::vector<int> const& getDetId2HitBitLUT() const { return mDetId2HitBitIndex; }
91+
8992
/// create a standalone ROOT file/tree with only the MCHeader branch
9093
static void extractFileFromKinematics(std::string_view kinefilename, std::string_view targetfilename);
9194

@@ -96,9 +99,10 @@ class MCEventHeader : public FairMCEventHeader
9699
// store a few global properties that this event
97100
// had in the current simulation (which can be used quick filtering/searching)
98101
MCEventStats mEventStats{};
102+
std::vector<int> mDetId2HitBitIndex;
99103
o2::utils::RootSerializableKeyValueStore mEventInfo;
100104

101-
ClassDefOverride(MCEventHeader, 3);
105+
ClassDefOverride(MCEventHeader, 4);
102106

103107
}; /** class MCEventHeader **/
104108

‎DataFormats/simulation/include/SimulationDataFormat/MCTrack.h‎

Lines changed: 27 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,8 @@ template <class _T>
4646
class MCTrackT
4747
{
4848
public:
49+
static constexpr int NHITBITS = 22; // do not modify this
50+
4951
/// Default constructor
5052
MCTrackT();
5153

@@ -159,24 +161,40 @@ class MCTrackT
159161
void SetSecondMotherTrackId(Int_t id) { mSecondMotherTrackId = id; }
160162
void SetFirstDaughterTrackId(Int_t id) { mFirstDaughterTrackId = id; }
161163
void SetLastDaughterTrackId(Int_t id) { mLastDaughterTrackId = id; }
164+
162165
// set bit indicating that this track
163-
// left a hit in detector with id iDet
164-
void setHit(Int_t iDet)
166+
// left a hit in detector that corresponds to bit iDetBit. This bit is found in
167+
// a lookup table.
168+
void setHit(Int_t iDetBit)
165169
{
166-
assert(0 <= iDet && iDet < o2::detectors::DetID::nDetectors);
170+
assert(0 <= iDetBit && iDetBit < o2::detectors::DetID::nDetectors);
167171
auto prop = ((PropEncoding)mProp);
168-
prop.hitmask |= 1 << iDet;
172+
prop.hitmask |= 1 << iDetBit;
169173
mProp = prop.i;
170174
}
171175

172-
// did detector iDet see this track?
173-
bool leftTrace(Int_t iDet) const { return (((PropEncoding)mProp).hitmask & (1 << iDet)) > 0; }
176+
bool leftTraceGivenBitField(int bit) const
177+
{
178+
return (((PropEncoding)mProp).hitmask & (1 << bit)) > 0;
179+
}
180+
181+
/// Returns of detector with DetID iDet has seen this track.
182+
/// Needs lookup table detIDToBit mapping detectorIDs to bitfields (which is persistified in MCEventHeaders).
183+
bool leftTrace(Int_t iDet, std::vector<int> const& detIDtoBit) const
184+
{
185+
auto bit = detIDtoBit[iDet];
186+
if (bit != -1) {
187+
return leftTraceGivenBitField(bit);
188+
}
189+
return false;
190+
}
191+
174192
// determine how many detectors "saw" this track
175193
int getNumDet() const
176194
{
177195
int count = 0;
178-
for (auto i = o2::detectors::DetID::First; i < o2::detectors::DetID::nDetectors; ++i) {
179-
if (leftTrace(i)) {
196+
for (auto i = 0; i < NHITBITS; ++i) {
197+
if (leftTraceGivenBitField(i)) {
180198
count++;
181199
}
182200
}
@@ -261,7 +279,7 @@ class MCTrackT
261279
struct {
262280
int storage : 1; // encoding whether to store this track to the output
263281
unsigned int process : 6; // encoding process that created this track (enough to store TMCProcess from ROOT)
264-
int hitmask : 22; // encoding hits per detector
282+
int hitmask : NHITBITS; // encoding hits per detector
265283
static_assert(o2::detectors::DetID::nDetectors <= 22); // ensure that all known detectors can be encoded here by a bit
266284
int reserved1 : 1; // bit reserved for possible future purposes
267285
int inhibited : 1; // whether tracking of this was inhibited

‎DataFormats/simulation/test/MCTrack.cxx‎

Lines changed: 11 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -24,10 +24,17 @@ using namespace o2;
2424
BOOST_AUTO_TEST_CASE(MCTrack_test)
2525
{
2626
MCTrack track;
27+
28+
// auxiliary lookup table needed to fetch and set hit properties
29+
std::vector<int> hitLUT(o2::detectors::DetID::nDetectors, -1);
30+
// in this test we have a single fictional detector 1, which we map to
31+
// the first bit
32+
hitLUT[1] = 0;
33+
2734
// check properties on default constructed object
2835
BOOST_CHECK(track.getStore() == false);
2936
for (auto i = o2::detectors::DetID::First; i < o2::detectors::DetID::nDetectors; ++i) {
30-
BOOST_CHECK(track.leftTrace(i) == false);
37+
BOOST_CHECK(track.leftTrace(i, hitLUT) == false);
3138
}
3239
BOOST_CHECK(track.getNumDet() == 0);
3340
BOOST_CHECK(track.hasHits() == false);
@@ -41,10 +48,10 @@ BOOST_AUTO_TEST_CASE(MCTrack_test)
4148
BOOST_CHECK(track.getStore() == true);
4249

4350
// set hit for first detector
44-
BOOST_CHECK(track.leftTrace(1) == false);
45-
track.setHit(1);
51+
BOOST_CHECK(track.leftTrace(1, hitLUT) == false);
52+
track.setHit(hitLUT[1]);
4653
BOOST_CHECK(track.hasHits() == true);
47-
BOOST_CHECK(track.leftTrace(1) == true);
54+
BOOST_CHECK(track.leftTrace(1, hitLUT) == true);
4855
BOOST_CHECK(track.getNumDet() == 1);
4956

5057
// check process encoding

‎Detectors/Base/include/DetectorsBase/Detector.h‎

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -216,6 +216,10 @@ class Detector : public FairDetector
216216
// which is required for media creation
217217
static void initFieldTrackingParams(int& mode, float& maxfield);
218218

219+
/// set the DetID to HitBitIndex mapping. Succeeds if not already set.
220+
static void setDetId2HitBitIndex(std::vector<int> const& v) { Detector::sDetId2HitBitIndex = v; }
221+
static std::vector<int> const& getDetId2HitBitIndex() { return Detector::sDetId2HitBitIndex; }
222+
219223
protected:
220224
Detector(const Detector& origin);
221225

@@ -233,6 +237,7 @@ class Detector : public FairDetector
233237

234238
static Float_t mDensityFactor; //! factor that is multiplied to all material densities (ONLY for
235239
// systematic studies)
240+
static std::vector<int> sDetId2HitBitIndex; //! global lookup table keeping mapping of DetID to index in hit bit field (used in MCTrack)
236241

237242
ClassDefOverride(Detector, 1); // Base class for ALICE Modules
238243
};
@@ -738,6 +743,7 @@ class DetImpl : public o2::base::Detector
738743
int mInitialized = false;
739744

740745
char* mHitCollectorBufferPtr = nullptr; //! pointer to hit (collector) buffer location (strictly internal)
746+
741747
ClassDefOverride(DetImpl, 0);
742748
};
743749
} // namespace base

‎Detectors/Base/src/Detector.cxx‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -30,6 +30,7 @@ using namespace o2::base;
3030
using namespace o2::detectors;
3131

3232
Float_t Detector::mDensityFactor = 1.0;
33+
std::vector<int> o2::base::Detector::sDetId2HitBitIndex{}; // initialize empty vector
3334

3435
Detector::Detector() : FairDetector(), mMapMaterial(), mMapMedium() {}
3536
Detector::Detector(const char* name, Bool_t Active)

‎Detectors/Base/src/Stack.cxx‎

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -638,14 +638,17 @@ void Stack::Print(Option_t* option) const
638638

639639
void Stack::addHit(int iDet)
640640
{
641+
// translate detector ID to bit id
642+
const auto bitID = o2::base::Detector::getDetId2HitBitIndex()[iDet];
643+
641644
if (mIndexOfCurrentTrack < mNumberOfPrimaryParticles) {
642645
auto& part = mTracks->at(mIndexOfCurrentTrack);
643-
part.setHit(iDet);
646+
part.setHit(bitID);
644647

645648
} else {
646649
Int_t iTrack = mTrackIDtoParticlesEntry[mIndexOfCurrentTrack];
647650
auto& part = mParticles[iTrack];
648-
part.setHit(iDet);
651+
part.setHit(bitID);
649652
}
650653
mCurrentParticle.SetBit(ParticleStatus::kHasHits, 1);
651654
mHitCounter++;
@@ -654,7 +657,9 @@ void Stack::addHit(int iDet, Int_t iTrack)
654657
{
655658
mHitCounter++;
656659
auto& part = mParticles[iTrack];
657-
part.setHit(iDet);
660+
// fetch the bit encoding for hits
661+
const auto bitID = o2::base::Detector::getDetId2HitBitIndex()[iDet];
662+
part.setHit(bitID);
658663
mCurrentParticle.SetBit(ParticleStatus::kHasHits, 1);
659664
}
660665

‎Steer/src/O2MCApplication.cxx‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -197,6 +197,7 @@ void O2MCApplicationBase::finishEventCommon()
197197

198198
auto header = static_cast<o2::dataformats::MCEventHeader*>(fMCEventHeader);
199199
header->getMCEventStats().setNSteps(mStepCounter);
200+
header->setDetId2HitBitLUT(o2::base::Detector::getDetId2HitBitIndex());
200201

201202
static_cast<o2::data::Stack*>(GetStack())->updateEventStats();
202203
}

‎macro/build_geometry.C‎

Lines changed: 36 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -179,103 +179,120 @@ void build_geometry(FairRunSim* run = nullptr)
179179
run->AddModule(new o2::passive::FrameStructure("FRAME", "Frame"));
180180
}
181181

182+
std::vector<int> detId2RunningId = std::vector<int>(o2::detectors::DetID::nDetectors, -1); // a mapping of detectorId to a dense runtime index
183+
// used for instance to set bits in the hit structure of MCTracks; -1 means that there is no bit associated
184+
185+
auto addReadoutDetector = [&detId2RunningId, &run](o2::base::Detector* detector) {
186+
static int runningid = 0; // this is static for constant lambda interfaces --> use fixed type and not auto in the lambda!
187+
run->AddModule(detector);
188+
if (detector->IsActive()) {
189+
auto detID = detector->GetDetId();
190+
detId2RunningId[detID] = runningid;
191+
LOG(info) << " DETID " << detID << " vs " << detector->GetDetId() << " mapped to hit bit index " << runningid;
192+
runningid++;
193+
}
194+
};
195+
182196
if (isActivated("TOF")) {
183197
// TOF
184-
run->AddModule(new o2::tof::Detector(isReadout("TOF")));
198+
addReadoutDetector(new o2::tof::Detector(isReadout("TOF")));
185199
}
186200

187201
if (isActivated("TRD")) {
188202
// TRD
189-
run->AddModule(new o2::trd::Detector(isReadout("TRD")));
203+
addReadoutDetector(new o2::trd::Detector(isReadout("TRD")));
190204
}
191205

192206
if (isActivated("TPC")) {
193207
// tpc
194-
run->AddModule(new o2::tpc::Detector(isReadout("TPC")));
208+
addReadoutDetector(new o2::tpc::Detector(isReadout("TPC")));
195209
}
196210
#ifdef ENABLE_UPGRADES
197211
if (isActivated("IT3")) {
198212
// IT3
199-
run->AddModule(new o2::its::Detector(isReadout("IT3"), "IT3"));
213+
addReadoutDetector(new o2::its::Detector(isReadout("IT3"), "IT3"));
200214
}
201215

202216
if (isActivated("TRK")) {
203217
// ALICE 3 TRK
204-
run->AddModule(new o2::trk::Detector(isReadout("TRK")));
218+
addReadoutDetector(new o2::trk::Detector(isReadout("TRK")));
205219
}
206220

207221
if (isActivated("FT3")) {
208222
// ALICE 3 FT3
209-
run->AddModule(new o2::ft3::Detector(isReadout("FT3")));
223+
addReadoutDetector(new o2::ft3::Detector(isReadout("FT3")));
210224
}
211225

212226
if (isActivated("FCT")) {
213227
// ALICE 3 FCT
214-
run->AddModule(new o2::fct::Detector(isReadout("FCT")));
228+
addReadoutDetector(new o2::fct::Detector(isReadout("FCT")));
215229
}
216230
#endif
217231

218232
if (isActivated("ITS")) {
219233
// its
220-
run->AddModule(new o2::its::Detector(isReadout("ITS")));
234+
addReadoutDetector(new o2::its::Detector(isReadout("ITS")));
221235
}
222236

223237
if (isActivated("MFT")) {
224238
// mft
225-
run->AddModule(new o2::mft::Detector(isReadout("MFT")));
239+
addReadoutDetector(new o2::mft::Detector(isReadout("MFT")));
226240
}
227241

228242
if (isActivated("MCH")) {
229243
// mch
230-
run->AddModule(new o2::mch::Detector(isReadout("MCH")));
244+
addReadoutDetector(new o2::mch::Detector(isReadout("MCH")));
231245
}
232246

233247
if (isActivated("MID")) {
234248
// mid
235-
run->AddModule(new o2::mid::Detector(isReadout("MID")));
249+
addReadoutDetector(new o2::mid::Detector(isReadout("MID")));
236250
}
237251

238252
if (isActivated("EMC")) {
239253
// emcal
240-
run->AddModule(new o2::emcal::Detector(isReadout("EMC")));
254+
addReadoutDetector(new o2::emcal::Detector(isReadout("EMC")));
241255
}
242256

243257
if (isActivated("PHS")) {
244258
// phos
245-
run->AddModule(new o2::phos::Detector(isReadout("PHS")));
259+
addReadoutDetector(new o2::phos::Detector(isReadout("PHS")));
246260
}
247261

248262
if (isActivated("CPV")) {
249263
// cpv
250-
run->AddModule(new o2::cpv::Detector(isReadout("CPV")));
264+
addReadoutDetector(new o2::cpv::Detector(isReadout("CPV")));
251265
}
252266

253267
if (isActivated("FT0")) {
254268
// FIT-T0
255-
run->AddModule(new o2::ft0::Detector(isReadout("FT0")));
269+
addReadoutDetector(new o2::ft0::Detector(isReadout("FT0")));
256270
}
257271

258272
if (isActivated("FV0")) {
259273
// FIT-V0
260-
run->AddModule(new o2::fv0::Detector(isReadout("FV0")));
274+
addReadoutDetector(new o2::fv0::Detector(isReadout("FV0")));
261275
}
262276

263277
if (isActivated("FDD")) {
264278
// FIT-FDD
265-
run->AddModule(new o2::fdd::Detector(isReadout("FDD")));
279+
addReadoutDetector(new o2::fdd::Detector(isReadout("FDD")));
266280
}
267281

268282
if (isActivated("HMP")) {
269283
// HMP
270-
run->AddModule(new o2::hmpid::Detector(isReadout("HMP")));
284+
addReadoutDetector(new o2::hmpid::Detector(isReadout("HMP")));
271285
}
272286

273287
if (isActivated("ZDC")) {
274288
// ZDC
275-
run->AddModule(new o2::zdc::Detector(isReadout("ZDC")));
289+
addReadoutDetector(new o2::zdc::Detector(isReadout("ZDC")));
276290
}
277291

278292
if (geomonly) {
279293
run->Init();
280294
}
295+
296+
// register the DetId2HitIndex lookup with the detector class by copying the vector
297+
o2::base::Detector::setDetId2HitBitIndex(detId2RunningId);
281298
}

0 commit comments

Comments
 (0)