Skip to content

Commit b3e5742

Browse files
shahor02shahor02
committed
Revert "QED events integrated in DigitizationContext"
This reverts commit b1d1bbf.
1 parent b1d1bbf commit b3e5742

5 files changed

Lines changed: 89 additions & 161 deletions

File tree

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

Lines changed: 8 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -26,19 +26,15 @@ namespace steer
2626
{
2727
// a structure describing EventPart
2828
// (an elementary constituent of a collision)
29-
30-
constexpr static int QEDSOURCEID = 99;
31-
3229
struct EventPart {
3330
EventPart() = default;
3431
EventPart(int s, int e) : sourceID(s), entryID(e) {}
3532
int sourceID = 0; // the ID of the source (0->backGround; > 1 signal source)
3633
// the sourceID should correspond to the chain ID
3734
int entryID = 0; // the event/entry ID inside the chain corresponding to sourceID
3835

39-
static bool isSignal(EventPart e) { return e.sourceID > 1 && e.sourceID != QEDSOURCEID; }
36+
static bool isSignal(EventPart e) { return e.sourceID > 1; }
4037
static bool isBackGround(EventPart e) { return !isSignal(e); }
41-
static bool isQED(EventPart e) { return e.sourceID == QEDSOURCEID; }
4238
ClassDefNV(EventPart, 1);
4339
};
4440

@@ -54,30 +50,24 @@ class DigitizationContext
5450
void setMaxNumberParts(int maxp) { mMaxPartNumber = maxp; }
5551
int getMaxNumberParts() const { return mMaxPartNumber; }
5652

57-
std::vector<o2::InteractionTimeRecord>& getEventRecords(bool withQED = false) { return withQED ? mEventRecordsWithQED : mEventRecords; }
58-
std::vector<std::vector<o2::steer::EventPart>>& getEventParts(bool withQED = false) { return withQED ? mEventPartsWithQED : mEventParts; }
59-
60-
const std::vector<o2::InteractionTimeRecord>& getEventRecords(bool withQED = false) const { return withQED ? mEventRecordsWithQED : mEventRecords; }
61-
const std::vector<std::vector<o2::steer::EventPart>>& getEventParts(bool withQED = false) const { return withQED ? mEventPartsWithQED : mEventParts; }
53+
std::vector<o2::InteractionTimeRecord>& getEventRecords() { return mEventRecords; }
54+
std::vector<std::vector<o2::steer::EventPart>>& getEventParts() { return mEventParts; }
6255

63-
bool isQEDProvided() const { return !mEventRecordsWithQED.empty(); }
56+
const std::vector<o2::InteractionTimeRecord>& getEventRecords() const { return mEventRecords; }
57+
const std::vector<std::vector<o2::steer::EventPart>>& getEventParts() const { return mEventParts; }
6458

6559
o2::BunchFilling& getBunchFilling() { return mBCFilling; }
6660
const o2::BunchFilling& getBunchFilling() const { return (const o2::BunchFilling&)mBCFilling; }
6761

6862
void setMuPerBC(float m) { mMuBC = m; }
6963
float getMuPerBC() const { return mMuBC; }
7064

71-
void printCollisionSummary(bool withQED = false) const;
65+
void printCollisionSummary() const;
7266

7367
// we need a method to fill the file names
7468
void setSimPrefixes(std::vector<std::string> const& p);
7569
std::vector<std::string> const& getSimPrefixes() const { return mSimPrefixes; }
7670

77-
/// add QED contributions to context; QEDprefix is prefix of QED production
78-
/// irecord is vector of QED interaction times (sampled externally)
79-
void fillQED(std::string_view QEDprefix, std::vector<o2::InteractionTimeRecord> const& irecord);
80-
8171
/// Common functions the setup input TChains for reading, given the state (prefixes) encapsulated
8272
/// by this context. The input vector needs to be empty otherwise nothing will be done.
8373
/// return boolean saying if input simchains was modified or not
@@ -117,17 +107,12 @@ class DigitizationContext
117107
// for each collision we record the constituents (which shall not exceed mMaxPartNumber)
118108
std::vector<std::vector<o2::steer::EventPart>> mEventParts;
119109

120-
// the collision records _with_ QED interleaved;
121-
std::vector<o2::InteractionTimeRecord> mEventRecordsWithQED;
122-
std::vector<std::vector<o2::steer::EventPart>> mEventPartsWithQED;
123-
124110
o2::BunchFilling mBCFilling; // patter of active BCs
125111

126-
std::vector<std::string> mSimPrefixes; // identifiers to the hit sim products; the key corresponds to the source ID of event record
127-
std::string mQEDSimPrefix; // prefix for QED production/contribution
112+
std::vector<std::string> mSimPrefixes; // identifiers to the hit sim products; the index corresponds to the source ID of event record
128113
mutable o2::parameters::GRPObject* mGRP = nullptr; //!
129114

130-
ClassDefNV(DigitizationContext, 3);
115+
ClassDefNV(DigitizationContext, 2);
131116
};
132117

133118
/// function reading the hits from a chain (previously initialized with initSimChains

‎DataFormats/simulation/src/DigitizationContext.cxx‎

Lines changed: 7 additions & 83 deletions
Original file line numberDiff line numberDiff line change
@@ -14,39 +14,25 @@
1414
#include <TChain.h>
1515
#include <TFile.h>
1616
#include <iostream>
17-
#include <numeric> // for iota
1817
#include <MathUtils/Cartesian3D.h>
1918

2019
using namespace o2::steer;
2120

22-
void DigitizationContext::printCollisionSummary(bool withQED) const
21+
void DigitizationContext::printCollisionSummary() const
2322
{
2423
std::cout << "Summary of DigitizationContext --\n";
2524
std::cout << "Parts per collision " << mMaxPartNumber << "\n";
2625
std::cout << "Collision parts taken from simulations specified by prefix:\n";
2726
for (int p = 0; p < mSimPrefixes.size(); ++p) {
2827
std::cout << "Part " << p << " : " << mSimPrefixes[p] << "\n";
2928
}
30-
if (withQED) {
31-
std::cout << "Number of Collisions " << mEventRecords.size() << "\n";
32-
std::cout << "Number of QED events " << mEventRecordsWithQED.size() - mEventRecords.size() << "\n";
33-
// loop over combined stuff
34-
for (int i = 0; i < mEventRecordsWithQED.size(); ++i) {
35-
std::cout << "Record " << i << " TIME " << mEventRecordsWithQED[i];
36-
for (auto& e : mEventPartsWithQED[i]) {
37-
std::cout << " (" << e.sourceID << " , " << e.entryID << ")";
38-
}
39-
std::cout << "\n";
40-
}
41-
} else {
42-
std::cout << "Number of Collisions " << mEventRecords.size() << "\n";
43-
for (int i = 0; i < mEventRecords.size(); ++i) {
44-
std::cout << "Collision " << i << " TIME " << mEventRecords[i];
45-
for (auto& e : mEventParts[i]) {
46-
std::cout << " (" << e.sourceID << " , " << e.entryID << ")";
47-
}
48-
std::cout << "\n";
29+
std::cout << "Number of Collisions " << mEventRecords.size() << "\n";
30+
for (int i = 0; i < mEventRecords.size(); ++i) {
31+
std::cout << "Collision " << i << " TIME " << mEventRecords[i];
32+
for (auto& e : mEventParts[i]) {
33+
std::cout << " (" << e.sourceID << " , " << e.entryID << ")";
4934
}
35+
std::cout << "\n";
5036
}
5137
}
5238

@@ -75,19 +61,6 @@ bool DigitizationContext::initSimChains(o2::detectors::DetID detid, std::vector<
7561
// add signal files
7662
simchains.back()->AddFile(o2::base::NameConf::getHitsFileName(detid, mSimPrefixes[source].data()).c_str());
7763
}
78-
79-
// QED part
80-
if (mEventRecordsWithQED.size() > 0) {
81-
if (mSimPrefixes.size() >= QEDSOURCEID) {
82-
LOG(FATAL) << "Too many signal chains; crashes with QED source ID";
83-
}
84-
85-
// it might be better to use an unordered_map for the simchains but this requires interface changes
86-
simchains.resize(QEDSOURCEID + 1, nullptr);
87-
simchains[QEDSOURCEID] = new TChain("o2sim");
88-
simchains[QEDSOURCEID]->AddFile(o2::base::NameConf::getHitsFileName(detid, mQEDSimPrefix).c_str());
89-
}
90-
9164
return true;
9265
}
9366

@@ -198,52 +171,3 @@ DigitizationContext const* DigitizationContext::loadFromFile(std::string_view fi
198171
file.GetObject("DigitizationContext", incontext);
199172
return incontext;
200173
}
201-
202-
void DigitizationContext::fillQED(std::string_view QEDprefix, std::vector<o2::InteractionTimeRecord> const& irecord)
203-
{
204-
mQEDSimPrefix = QEDprefix;
205-
206-
std::vector<std::vector<o2::steer::EventPart>> qedEventParts;
207-
208-
// we need to fill the QED parts (using a simple round robin scheme)
209-
auto qedKinematicsName = o2::base::NameConf::getMCKinematicsFileName(mQEDSimPrefix);
210-
// find out how many events are stored
211-
TFile f(qedKinematicsName.c_str(), "OPEN");
212-
auto t = (TTree*)f.Get("o2sim");
213-
if (!t) {
214-
LOG(ERROR) << "No QED kinematics found";
215-
return;
216-
}
217-
auto numberQEDevents = t->GetEntries();
218-
int eventID = 0;
219-
for (auto& tmp : irecord) {
220-
std::vector<EventPart> qedpart;
221-
qedpart.emplace_back(QEDSOURCEID, eventID++);
222-
qedEventParts.push_back(qedpart);
223-
if (eventID == numberQEDevents) {
224-
eventID = 0;
225-
}
226-
}
227-
228-
// we need to do the interleaved event records for detectors consuming both
229-
// normal and QED events
230-
// --> merge both; sort first according to times and sort second one according to same order
231-
auto combinedrecords = mEventRecords;
232-
combinedrecords.insert(combinedrecords.end(), irecord.begin(), irecord.end());
233-
auto combinedparts = mEventParts;
234-
combinedparts.insert(combinedparts.end(), qedEventParts.begin(), qedEventParts.end());
235-
236-
// get sorted index vector based on event records
237-
std::vector<size_t> idx(combinedrecords.size());
238-
std::iota(idx.begin(), idx.end(), 0);
239-
240-
std::stable_sort(idx.begin(), idx.end(),
241-
[&combinedrecords](size_t i1, size_t i2) { return combinedrecords[i1] < combinedrecords[i2]; });
242-
243-
mEventRecordsWithQED.clear();
244-
mEventPartsWithQED.clear();
245-
for (int i = 0; i < idx.size(); ++i) {
246-
mEventRecordsWithQED.push_back(combinedrecords[idx[i]]);
247-
mEventPartsWithQED.push_back(combinedparts[idx[i]]);
248-
}
249-
}

‎Steer/DigitizerWorkflow/src/ITSMFTDigitizerSpec.cxx‎

Lines changed: 72 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,13 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
4848
using BaseDPLDigitizer::init;
4949
void initDigitizerTask(framework::InitContext& ic) override
5050
{
51+
// init optional QED chain
52+
auto qedfilename = ic.options().get<std::string>("simFileQED");
53+
if (qedfilename.size() > 0) {
54+
mQEDChain.AddFile(qedfilename.c_str());
55+
LOG(INFO) << "Attach QED Tree: " << mQEDChain.GetEntries();
56+
}
57+
5158
setDigitizationOptions(); // set options provided via configKeyValues mechanism
5259
auto& digipar = mDigitizer.getParams();
5360

@@ -81,11 +88,8 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
8188
// read collision context from input
8289
auto context = pc.inputs().get<o2::steer::DigitizationContext*>("collisioncontext");
8390
context->initSimChains(mID, mSimChains);
84-
const bool withQED = context->isQEDProvided();
85-
auto& timesview = context->getEventRecords(withQED);
86-
LOG(INFO) << "GOT " << timesview.size() << " COLLISSION TIMES";
87-
LOG(INFO) << "SIMCHAINS " << mSimChains.size();
88-
91+
auto& timesview = context->getEventRecords();
92+
LOG(DEBUG) << "GOT " << timesview.size() << " COLLISSION TIMES";
8993
// if there is nothing to do ... return
9094
if (timesview.size() == 0) {
9195
return;
@@ -98,11 +102,18 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
98102
mDigitizer.setROFRecords(&mROFRecords);
99103
mDigitizer.setMCLabels(&mLabels);
100104

101-
auto& eventParts = context->getEventParts(withQED);
105+
// attach optional QED digits branch
106+
setupQEDChain();
107+
108+
auto& eventParts = context->getEventParts();
102109
// loop over all composite collisions given from context (aka loop over all the interaction records)
103110
for (int collID = 0; collID < timesview.size(); ++collID) {
104111
const auto& irt = timesview[collID];
105112

113+
if (mQEDChain.GetEntries()) { // QED must be processed before other inputs since done in small time steps
114+
processQED(irt);
115+
}
116+
106117
mDigitizer.setEventTime(irt);
107118
mDigitizer.resetEventROFrames(); // to estimate min/max ROF for this collID
108119
// for each collision, loop over the constituents event and source IDs
@@ -113,15 +124,19 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
113124
mHits.clear();
114125
context->retrieveHits(mSimChains, o2::detectors::SimTraits::DETECTORBRANCHNAMES[mID][0].c_str(), part.sourceID, part.entryID, &mHits);
115126

116-
if (mHits.size() > 0) {
117-
LOG(INFO) << "For collision " << collID << " eventID " << part.entryID
118-
<< " found " << mHits.size() << " hits ";
119-
mDigitizer.process(&mHits, part.entryID, part.sourceID); // call actual digitization procedure
120-
}
127+
LOG(INFO) << "For collision " << collID << " eventID " << part.entryID
128+
<< " found " << mHits.size() << " hits ";
129+
130+
mDigitizer.process(&mHits, part.entryID, part.sourceID); // call actual digitization procedure
121131
}
122132
mMC2ROFRecordsAccum.emplace_back(collID, -1, mDigitizer.getEventROFrameMin(), mDigitizer.getEventROFrameMax());
123133
accumulate();
124134
}
135+
// finish digitization ... stream any remaining digits/labels
136+
if (mQEDChain.GetEntries()) { // fill last slots from QED input
137+
processQED(mDigitizer.getEndTimeOfROFMax());
138+
}
139+
125140
mDigitizer.fillOutputContainer();
126141
accumulate();
127142

@@ -147,6 +162,42 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
147162
protected:
148163
ITSMFTDPLDigitizerTask(bool mctruth = true) : BaseDPLDigitizer(InitServices::FIELD | InitServices::GEOM), mWithMCTruth(mctruth) {}
149164

165+
void processQED(const o2::InteractionTimeRecord& irt)
166+
{
167+
168+
auto tQEDNext = mLastQEDTime.getTimeNS() + mQEDEntryTimeBinNS; // timeslice to retrieve
169+
std::string detStr = mID.getName();
170+
auto br = mQEDChain.GetBranch((detStr + "Hit").c_str());
171+
while (tQEDNext < irt.getTimeNS()) {
172+
mLastQEDTime.setFromNS(tQEDNext); // time used for current QED slot
173+
tQEDNext += mQEDEntryTimeBinNS; // prepare time for next QED slot
174+
if (++mLastQEDEntry >= mQEDChain.GetEntries()) {
175+
mLastQEDEntry = 0; // wrapp if needed
176+
}
177+
br->GetEntry(mLastQEDEntry);
178+
mDigitizer.setEventTime(mLastQEDTime);
179+
mDigitizer.process(&mHits, mLastQEDEntry, QEDSourceID);
180+
//
181+
}
182+
}
183+
184+
void setupQEDChain()
185+
{
186+
if (!mQEDChain.GetEntries()) {
187+
return;
188+
}
189+
mLastQEDTime.setFromNS(0.);
190+
std::string detStr = mID.getName();
191+
auto qedBranch = mQEDChain.GetBranch((detStr + "Hit").c_str());
192+
assert(qedBranch != nullptr);
193+
assert(mQEDEntryTimeBinNS >= 1.0);
194+
assert(QEDSourceID < o2::MCCompLabel::maxSourceID());
195+
mLastQEDTimeNS = -mQEDEntryTimeBinNS / 2; // time will be assigned to the middle of the bin
196+
qedBranch->SetAddress(&mHitsP);
197+
LOG(INFO) << "Attaching QED ITS hits as sourceID=" << int(QEDSourceID) << ", entry integrates "
198+
<< mQEDEntryTimeBinNS << " ns";
199+
}
200+
150201
void accumulate()
151202
{
152203
// accumulate result of single event processing, called after processing every event supplied
@@ -204,9 +255,16 @@ class ITSMFTDPLDigitizerTask : BaseDPLDigitizer
204255
o2::dataformats::MCTruthContainer<o2::MCCompLabel> mLabelsAccum;
205256
std::vector<o2::itsmft::MC2ROFRecord> mMC2ROFRecordsAccum;
206257
std::vector<TChain*> mSimChains;
258+
TChain mQEDChain = {"o2sim"};
207259

260+
double mQEDEntryTimeBinNS = 1000; // time-coverage of single QED tree entry in ns (TODO: make it settable)
261+
o2::InteractionTimeRecord mLastQEDTime;
262+
double mLastQEDTimeNS = 0; // time assingned to last QED entry
263+
int mLastQEDEntry = -1; // last used QED entry
208264
int mFixMC2ROF = 0; // 1st entry in mc2rofRecordsAccum to be fixed for ROFRecordID
209265
o2::parameters::GRPObject::ROMode mROMode = o2::parameters::GRPObject::PRESENT; // readout mode
266+
267+
const int QEDSourceID = 99; // unique source ID for the QED (TODO: move it as a const to general class?)
210268
};
211269

212270
//_______________________________________________
@@ -322,6 +380,7 @@ DataProcessorSpec getITSDigitizerSpec(int channel, bool mctruth)
322380
makeOutChannels(detOrig, mctruth),
323381
AlgorithmSpec{adaptFromTask<ITSDPLDigitizerTask>(mctruth)},
324382
Options{
383+
{"simFileQED", VariantType::String, "", {"Sim (QED) input filename"}},
325384
// { "configKeyValues", VariantType::String, "", { parHelper.str().c_str() } }
326385
}};
327386
}
@@ -341,7 +400,8 @@ DataProcessorSpec getMFTDigitizerSpec(int channel, bool mctruth)
341400
static_cast<SubSpecificationType>(channel), Lifetime::Timeframe}},
342401
makeOutChannels(detOrig, mctruth),
343402
AlgorithmSpec{adaptFromTask<MFTDPLDigitizerTask>(mctruth)},
344-
Options{}};
403+
Options{
404+
{"simFileQED", VariantType::String, "", {"Sim (QED) input filename"}}}};
345405
}
346406

347407
} // end namespace itsmft

0 commit comments

Comments
 (0)