8 #include "CLHEP/Random/RandGaussZiggurat.h"
35 if (m_runVoltage < 2.3 || m_runVoltage > 3.2){
36 ATH_MSG_WARNING(
"STGC run voltage must be within fit domain of 2.3 kV to 3.2 kV");
37 return StatusCode::FAILURE;
44 return StatusCode::SUCCESS;
49 int bunchInteger = (digitTime > 0) ?
static_cast<int>(std::abs(digitTime / 25.0))
50 :
static_cast<int>(std::abs(digitTime / 25.0)) + 1;
51 bctag = (bctag | bunchInteger);
52 if (digitTime < 0) bctag = ~bctag;
60 float elecThreshold = 0.0;
61 if (!thresholdData.
getThreshold(channelID, elecThreshold)) {
62 THROW_EXCEPTION(
"Cannot retrieve VMM threshold from conditions database!");
66 THROW_EXCEPTION(
"Cannot convert VMM charge threshold via conditions data!");
89 const double angular_tolerance = 1
e-3;
92 if (viewer.size() == 0) {
96 std::array<sTgcSimDigitVec , 3> simDigitsByChType{};
104 double eventTime = hit.eventTime();
105 earliestEventTime =
std::min(earliestEventTime, eventTime);
113 const Amg::Vector3D locHitPos = xAOD::toEigen(hit->localPosition());
114 const Amg::Vector3D locHitDir = xAOD::toEigen(hit->localDirection());
116 const double hitKineticEnergy = hit->kineticEnergy();
118 ATH_MSG_DEBUG(
"Skip electron hit with kinetic energy " << hitKineticEnergy
125 const double theta = std::acos(locHitDir.z());
128 if (std::abs(theta - ninetyDegrees) < angular_tolerance) {
129 ATH_MSG_VERBOSE(
"Skipping hit nearly perpendicular to Z-axis (angle: " << theta <<
")");
132 if (std::abs(locHitDir.z()) < 0.00001) {
137 ATH_MSG_DEBUG(
"Updated hit global time to include off set of " << eventTime <<
" ns from OOT bunch.");
147 bool acceptHit =
true;
156 double globalHitTime = hit->globalTime() + eventTime;
158 double bunchTime = globalHitTime - tofCorrection;
163 <<
" globalHitPosition " <<
Amg::toString(globalHitPos) <<
" hit: r " << globalHitPos.perp() <<
" z " << globalHitPos.z()
164 <<
" mclink " << particleLink <<
"Kinetic energy " << hitKineticEnergy);
165 ATH_MSG_VERBOSE(
"Total hits passed to digitize: " << hitsToDigit.size());
167 if(digitizedHits.empty()) {
170 ATH_MSG_VERBOSE(
"Hit produced " << digitizedHits.size() <<
" digits." );
172 for (std::unique_ptr<sTgcDigit>&
digit : digitizedHits) {
182 double digitTime =
digit->time();
183 int digitChType = idHelper.channelType(digitId);
185 if(digitChType == ReadoutChannelType::Strip) {
195 double digitCharge =
digit->charge();
197 int eventId = hit.eventId();
198 bool isDead{
false},
isPileup{eventId != 0};
199 ATH_MSG_VERBOSE(
"Hit is from the main signal subevent if eventId is zero, eventId = "
200 << eventId <<
" newDigit time: " << digitTime);
201 auto newDigitPtr = std::make_unique<sTgcDigit>(digitId, digitBCTag, digitTime, digitCharge, isDead,
isPileup);
202 if (digitChType == ReadoutChannelType::Strip) {
204 <<
" BC tag = " << newDigitPtr->bcTag()
205 <<
" digitTime = " << newDigitPtr->time()
206 <<
" charge = " << newDigitPtr->charge());
208 simDigitsByChType[digitChType].emplace_back(std::move(hit), std::move(newDigitPtr));
211 if (!simDigitsByChType[ReadoutChannelType::Strip].empty()) {
214 if (!simDigitsByChType[ReadoutChannelType::Pad].empty()) {
217 if (!simDigitsByChType[ReadoutChannelType::Wire].empty()) {
220 }
while (viewer.next());
224 return StatusCode::SUCCESS;
230 const double vmmDeadTime,
231 const bool isNeighbourOn,
236 if (digitsInChamber.empty()) {
237 ATH_MSG_WARNING(
"Failed to obtain the digitized hits for VMM Simulation" );
238 return StatusCode::SUCCESS;
242 digitsInChamber, isNeighbourOn);
244 if (mergedDigits.empty()) {
245 return StatusCode::SUCCESS;
251 bool acceptDigit{
true};
252 float chargeAfterSmearing = merged.getDigit().charge();
262 if (idHelper.channelType(merged.identify()) == ReadoutChannelType::Strip &&
263 chargeAfterSmearing < 0.001) {
266 std::unique_ptr<sTgcDigit> finalDigit = std::make_unique<sTgcDigit>(std::move(merged.getDigit()));
271 " BC tag = " << finalDigit->
bcTag()<<
272 " digitTime = " << finalDigit->
time() <<
273 " charge = " << finalDigit->
charge());
279 double globalHitTime = sdoHit->
globalTime() + merged.getSimHit().eventTime();
284 outColl->
push_back(std::move(finalDigit));
286 return StatusCode::SUCCESS;
291 const double vmmDeadTime,
293 const bool isNeighbourOn)
const {
298 std::stable_sort(unmergedDigits.begin(), unmergedDigits.end(),
300 const int layA = idHelper.gasGap(a.identify());
301 const int layB = idHelper.gasGap(b.identify());
302 if (layA != layB) return layA < layB;
303 const int chA = idHelper.channel(a.identify());
304 const int chB = idHelper.channel(b.identify());
305 if (chA != chB) return chA < chB;
306 return a.time() < b.time();
311 premerged.reserve(unmergedDigits.size());
312 savedDigits.reserve(premerged.capacity());
315 if (!isNeighbourOn || savedDigits.empty())
return false;
316 if (savedDigits.back().identify() == candidate.identify() &&
317 std::abs(savedDigits.back().time() - candidate.time()) < vmmDeadTime) {
321 const Identifier digitId = candidate.identify();
322 const int channel = idHelper.channel(digitId);
323 const int maxChannel = detMgr->getsTgcReadoutElement(digitId)->numChannels(digitId);
326 if (neighbour ==
channel)
continue;
327 const Identifier neighbourId = idHelper.channelID(digitId,
328 idHelper.multilayer(digitId),
329 idHelper.gasGap(digitId),
330 idHelper.channelType(digitId), neighbour);
334 return known.identify() == neighbourId &&
335 known.getDigit().charge() > threshold &&
336 std::abs(known.time() - candidate.time()) < m_hitTimeMergeThreshold;
337 }) != savedDigits.end())
return true;
351 sTgcDigit& digit1{(*merge_me).getDigit()};
352 double totalCharge = digit1.
charge();
353 double weightedTime = digit1.time();
356 for ( ; merge_with!= unmergedDigits.end(); ++merge_with) {
358 if ((*merge_with).identify() != (*merge_me).identify()) {
361 const sTgcDigit& mergeDigit{(*merge_with).getDigit()};
367 weightedTime = (weightedTime * totalCharge + mergeDigit.time() * mergeDigit.charge())
368 / (totalCharge + mergeDigit.charge());
370 totalCharge += mergeDigit.
charge();
372 digit1.set_charge(totalCharge);
373 digit1.set_time(weightedTime);
375 if (!savedDigits.empty() &&
376 savedDigits.back().identify() == digit1.identify() &&
377 std::abs(savedDigits.back().time() - digit1.time()) <= vmmDeadTime)
continue;
378 if (digit1.charge() >
threshold || passNeigbourLogic(mergedHit)){
379 savedDigits.emplace_back(std::move(mergedHit));
380 }
else if (isNeighbourOn) {
381 premerged.emplace_back(std::move(mergedHit));
384 std::copy_if(std::make_move_iterator(premerged.begin()),
385 std::make_move_iterator(premerged.end()),
386 std::back_inserter(savedDigits), passNeigbourLogic);
387 if(savedDigits.empty() && !unmergedDigits.empty()) {