Digitization functionality shared with RPC_PileUpTool.
Use special jitter consant for BIS & BIL chambers.
486 {
487 ATHRNG::RNGWrapper* rngWrapper =
m_rndmSvc->getEngine(
this);
489 CLHEP::HepRandomEngine* rndmEngine = rngWrapper->
getEngine(ctx);
490
491 const MuonGM::MuonDetectorManager* detMgr{nullptr};
493
494
495 std::unique_ptr<RPCSimHitCollection> inputSimHitColl{std::make_unique<RPCSimHitCollection>("RPC_Hits")};
496
497
498
499
501
502
505 return StatusCode::FAILURE;
506 }
507
508 struct SimDataContent {
510 std::vector<MuonSimData::Deposit> deposits;
512 double simTime{0.};
513 };
514
515 while (
m_thpcRPC->nextDetectorElement(i, e)) {
516
517
518 std::map<Identifier, SimDataContent> channelSimDataMap;
519
520
521 while (i != e) {
523
524 TimedHitPtr<RPCSimHit> phit(*i++);
525
526
527 const RPCSimHit&
hit(*phit);
528
529 const int idHit =
hit.RPCid();
530
531 const double globalHitTime{
hitTime(phit)};
532
533 const double G4Time{
hit.globalTime()};
534
535 const double bunchTime{globalHitTime -
hit.globalTime()};
536
537 ATH_MSG_DEBUG(
"Global time " << globalHitTime <<
" G4 time " << G4Time <<
" Bunch time " << bunchTime);
538
540 ATH_MSG_VERBOSE(
"Validation: globalHitTime, G4Time, BCtime = " << globalHitTime <<
" " << G4Time <<
" " << bunchTime);
541 inputSimHitColl->Emplace(idHit, globalHitTime,
hit.localPosition(),
543 hit.postLocalPosition(),
544 hit.energyDeposit(),
hit.stepLength(),
hit.particleEncoding(),
hit.kineticEnergy());
545 }
546
547
555
557
558
560 const Identifier elementID =
m_idHelper->elementID(stationName,stationEta,stationPhi,doubletR,
isValid);
562 ATH_MSG_WARNING(
"Failed to construct the element ID from "<<stationName
563 <<", stationEta: "<<stationEta<<", stationPhi: "<<stationPhi<<", doubletR: "<<doubletR);
564 continue;
565 }
566
568 << " stationName " << stationName << " stationEta " << stationEta << " stationPhi " << stationPhi << " doubletR "
569 << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi << " gasGap " << gasGap);
570 const Identifier detElId{
m_idHelper->channelID(elementID, doubletZ, doubletPhi, 1,0, 1,
isValid)};
572 continue;
573 }
578 }
579
580
581 bool isValidEta{false}, isValidPhi{false};
582 const Identifier idpaneleta =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 0, 1, isValidEta);
583 const Identifier idpanelphi =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 1, 1, isValidPhi);
584 if (!isValidEta || !isValidPhi) {
586 << " stationName " << stationName << " stationEta " << stationEta << " stationPhi " << stationPhi
587 << " doubletR " << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi << " gasGap "
588 << gasGap);
589 continue;
590 }
591
592
596 const double corrtimejitter = tmp_CorrJitter > 0.01 ?
597 CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0., tmp_CorrJitter) : 0.;
598
599
600
601
602
606
607 std::array<int, 3> pcseta =
physicalClusterSize(ctx, reEle, idpaneleta, gapCentre, rndmEngine);
608 ATH_MSG_VERBOSE(
"Simulated cluster on eta panel: size/first/last= " << pcseta[0] <<
"/" << pcseta[1] <<
"/" << pcseta[2]);
609 std::array<int, 3> pcsphi =
physicalClusterSize(ctx, reEle, idpanelphi, gapCentre, rndmEngine);
610 ATH_MSG_VERBOSE(
"Simulated cluster on phi panel: size/first/last= " << pcsphi[0] <<
"/" << pcsphi[1] <<
"/" << pcsphi[2]);
611
612
613
614
615 const Identifier atlasRpcIdeta =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 0, pcseta[1], isValidEta);
616 const Identifier atlasRpcIdphi =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 1, pcsphi[1], isValidPhi);
617
619 const auto [etaStripOn, phiStripOn] =
detectionEfficiency(ctx, idpaneleta, idpanelphi, rndmEngine, particleLink);
620 ATH_MSG_DEBUG(
"SetPhiOn " << phiStripOn <<
" SetEtaOn " << etaStripOn);
621
622 for (bool imeasphi : {false, true}) {
623 if (!imeasphi && (!etaStripOn || !isValidEta)) continue;
624 if (imeasphi && (!phiStripOn || !isValidPhi)) continue;
625
626
627
628 const Identifier& atlasId = !imeasphi ? atlasRpcIdeta : atlasRpcIdphi;
629 std::array<int, 3>
pcs{!imeasphi ? pcseta : pcsphi};
630
631 ATH_MSG_DEBUG(
"SetOn: stationName " << stationName <<
" stationEta " << stationEta <<
" stationPhi " << stationPhi
632 << " doubletR " << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi
633 << " gasGap " << gasGap << " measphi " << imeasphi);
634
635
637 if (pcs[2] < 0){
638 continue;
639 }
640
641 ATH_MSG_DEBUG(
"Simulated cluster1: size/first/last= " << pcs[0] <<
"/" << pcs[1] <<
"/" << pcs[2]);
642
643
646
648 <<
" hit "<<
m_idHelper->print_to_string(atlasId)
656
657
659
660 double tns = G4Time + proptime + corrtimejitter;
661 ATH_MSG_VERBOSE(
"TOF+propagation time " << tns <<
" /s where proptime " << proptime <<
"/s");
662
663 double time = tns + bunchTime;
665
666
668
669 double*
b =
reinterpret_cast<double*
>(&packedMCword);
670
672
673
674
676
677
678
680 if (std::abs(
hit.particleEncoding()) == 13 ||
hit.particleEncoding() == 0) {
681 if (channelSimDataMap.find(atlasId) == channelSimDataMap.end()) {
682 SimDataContent&
content = channelSimDataMap[atlasId];
684 content.deposits.push_back(deposit);
689 }
690 }
691 }
692
693
694
695
696
697
698
700 for (
int clus = pcs[1]; clus <=
pcs[2]; ++clus) {
701 Identifier newId =
m_idHelper->channelID(stationName, stationEta, stationPhi, doubletR, doubletZ,
702 doubletPhi, gasGap, imeasphi, clus,
isValid);
704 ATH_MSG_WARNING(__FILE__<<
":"<<__LINE__<<
"Channel "<< stationName<<
" "<<stationEta<<
" "<<stationPhi<<
" "<< doubletR<<
" "<<doubletZ
705 <<" "<< doubletPhi<<" "<< gasGap <<" "<< imeasphi<<" "<< clus<<" is invalid");
706 continue;
707 }
708
711 ATH_MSG_WARNING(
"Temporary skipping creation of RPC digit for stationName="
712 << stationName << ", eta=" << stationEta << ", phi=" << stationPhi << ", doubletR=" << doubletR
713 << ", doubletZ=" << doubletZ << ", doubletPhi=" << doubletPhi << ", gasGap=" << gasGap
714 << ", measuresPhi=" << imeasphi << ", strip=" << clus << ", cf. ATLASRECTS-6124");
715 return StatusCode::SUCCESS;
716 } else {
719 return StatusCode::FAILURE;
720 }
721 }
722
726
728 std::vector<MuonSimData::Deposit> newdeps;
729 newdeps.push_back(deposit);
730 m_sdo_tmp_map.insert(std::map<Identifier, std::vector<MuonSimData::Deposit>>::value_type(newId, newdeps));
731 } else {
733 }
734 }
735 }
736 }
737
739 for (auto it = channelSimDataMap.begin(); it != channelSimDataMap.end(); ++it) {
740 MuonSimData
simData(
it->second.deposits, 0);
743 auto insertResult = sdoContainer->insert(std::make_pair(
it->first,
simData));
744 if (!insertResult.second)
745 ATH_MSG_WARNING(
"Attention: this sdo is not recorded, since the identifier already exists in the sdoContainer map");
746 }
747 }
748
749 }
750
752
753 std::map<Identifier, std::vector<MuonSimData::Deposit>>::iterator map_iter =
m_sdo_tmp_map.begin();
755
757
758 const Identifier theId = (*map_iter).first;
760
761 const std::vector<MuonSimData::Deposit> theDeps = (*map_iter).second;
762
763
765 std::multimap<double, MuonSimData::Deposit>
times;
766
767
768 for (
unsigned int k = 0;
k < theDeps.size();
k++) {
769 double time = theDeps[
k].second.secondEntry();
770 times.insert(std::multimap<double, MuonSimData::Deposit>::value_type(time, theDeps[k]));
771 }
772
773
774
775 IdContext rpcContext =
m_idHelper->module_context();
776
777 std::multimap<double, MuonSimData::Deposit>::iterator map_dep_iter =
times.begin();
778
779
780 double last_time = -10000;
781 for (; map_dep_iter !=
times.end(); ++map_dep_iter) {
782 double currTime = (*map_dep_iter).first;
784
786
787
788 if (sdoContainer->find(theId) != sdoContainer->end())
789 {
790 std::map<Identifier, MuonSimData>::const_iterator
it = sdoContainer->find(theId);
791 std::vector<MuonSimData::Deposit> deps = ((*it).second).getdeposits();
792 deps.push_back((*map_dep_iter).second);
793 } else
794 {
795 std::vector<MuonSimData::Deposit> deposits;
796 deposits.push_back((*map_dep_iter).second);
797 std::pair<std::map<Identifier, MuonSimData>::iterator, bool> insertResult =
798 sdoContainer->insert(std::make_pair(theId, MuonSimData(deposits, 0)));
799 if (!insertResult.second)
801 "Attention TEMP: this sdo is not recorded, since the identifier already exists in the sdoContainer map");
802 }
803 }
804
805 if (std::abs(currTime - last_time) > (
m_deadTime)) {
806 ATH_MSG_DEBUG(
"deposit with time " << currTime <<
" is distant enough from previous (if any) hit on teh same strip");
807 last_time = (*map_dep_iter).first;
808
809
810 double uncorrjitter = 0;
813 if (tmp_UncorrJitter > 0.01) uncorrjitter = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0., tmp_UncorrJitter);
814
815
816
817
821
823 double newDigit_time = currTime + uncorrjitter +
m_rpc_time_shift - tp - propTimeFromStripCenter;
824
825 double digi_ToT = -1.;
827
828 ATH_MSG_VERBOSE(
"last_time=currTime " << last_time <<
" jitter " << uncorrjitter <<
" TOFcorrection " << tp <<
" shift "
830
831
832 bool outsideDigitizationWindow =
outsideWindow(newDigit_time);
833 if (outsideDigitizationWindow) {
834 ATH_MSG_VERBOSE(
"hit outside digitization window - do not produce digits");
838
839 continue;
840 }
841
842
843 last_time = (*map_dep_iter).first;
844
845 std::unique_ptr<RpcDigit> newDigit = std::make_unique<RpcDigit>(theId, newDigit_time, digi_ToT, false);
846
847 Identifier elemId =
m_idHelper->elementID(theId);
848 RpcDigitCollection* digitCollection = nullptr;
849
850 IdentifierHash coll_hash;
851 if (
m_idHelper->get_hash(elemId, coll_hash, &rpcContext)) {
852 ATH_MSG_ERROR(
"Unable to get RPC hash id from RPC Digit collection "
853 <<
"context begin_index = " << rpcContext.
begin_index()
854 <<
" context end_index = " << rpcContext.
end_index() <<
" the identifier is \n"<<elemId);
855 }
856
857
859
860
861 if (coll_hash >= collections.size()) {
862 collections.resize (coll_hash+1);
863 }
864 digitCollection = collections[coll_hash].
get();
865 if (!digitCollection) {
866 collections[coll_hash] = std::make_unique<RpcDigitCollection>(elemId, coll_hash);
867 digitCollection = collections[coll_hash].
get();
868 }
869 digitCollection->
push_back(std::move(newDigit));
870
872
873 if (sdoContainer->find(theId) != sdoContainer->end()) {
874 std::map<Identifier, MuonSimData>::const_iterator
it = sdoContainer->find(theId);
875 std::vector<MuonSimData::Deposit> deps = ((*it).second).getdeposits();
876 deps.push_back((*map_dep_iter).second);
877 } else {
878 std::vector<MuonSimData::Deposit> deposits;
879 deposits.push_back((*map_dep_iter).second);
880 std::pair<std::map<Identifier, MuonSimData>::iterator, bool> insertResult =
881 sdoContainer->insert(std::make_pair(theId, MuonSimData(deposits, 0)));
882 if (!insertResult.second)
884 "Attention: this sdo is not recorded, since teh identifier already exists in the sdoContainer map");
885 }
886 }
887
888 } else
889 ATH_MSG_DEBUG(
"discarding digit due to dead time: " << (*map_dep_iter).first <<
" " << last_time);
890 }
891
892 }
893
894
899 }
900
901 return StatusCode::SUCCESS;
902}
float hitTime(const AFP_SIDSimHit &hit)
#define ATH_CHECK
Evaluate an expression and check for errors.
bool isValid() const
Test to see if the link can be dereferenced.
void setSeed(const std::string &algName, const EventContext &ctx)
Set the random seed using a string (e.g.
CLHEP::HepRandomEngine * getEngine(const EventContext &ctx) const
Retrieve the random engine corresponding to the provided EventContext.
const T * get(size_type n) const
Access an element, as an rvalue.
value_type push_back(value_type pElem)
Add an element to the end of the collection.
static HepMcParticleLink getRedirectedLink(const HepMcParticleLink &particleLink, uint32_t eventIndex, const EventContext &ctx)
Return a HepMcParticleLink pointing at the same particle, but in a different GenEvent.
size_type begin_index() const
size_type end_index() const
virtual const Amg::Transform3D & transform() const override
Return local to global transform.
const RpcReadoutElement * getRpcReadoutElement(const Identifier &id) const
access via extended identifier (requires unpacking)
bool rotatedRpcModule() const
Amg::Vector3D stripPos(const Identifier &id) const
Amg::Transform3D localToGlobalTransf(const Identifier &id) const
std::pair< HepMcParticleLink, MuonMCData > Deposit
StatusCode record(std::unique_ptr< T > data)
Record a const object to the store.
TimedVector::const_iterator const_iterator
std::optional< double > intersect(const AmgVector(N)&posA, const AmgVector(N)&dirA, const AmgVector(N)&posB, const AmgVector(N)&dirB)
Calculates the point B' along the line B that's closest to a second line A.
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
Eigen::Matrix< double, 3, 1 > Vector3D
time(flags, cells_name, *args, **kw)
bool ignoreTruthLink(const T &p, bool vetoPileUp)
Helper function for SDO creation in PileUpTools.
constexpr uint8_t stationPhi
station Phi 1 to 8