Digitization functionality shared with RPC_PileUpTool.
Use special jitter consant for BIS & BIL chambers.
482 {
483 ATHRNG::RNGWrapper* rngWrapper =
m_rndmSvc->getEngine(
this);
485 CLHEP::HepRandomEngine* rndmEngine = rngWrapper->
getEngine(ctx);
486
487 const MuonGM::MuonDetectorManager* detMgr{nullptr};
489
490
491 std::unique_ptr<RPCSimHitCollection> inputSimHitColl{std::make_unique<RPCSimHitCollection>("RPC_Hits")};
492
493
494
495
497
498
501 return StatusCode::FAILURE;
502 }
503
504 struct SimDataContent {
506 std::vector<MuonSimData::Deposit> deposits;
508 double simTime{0.};
509 };
510
511 while (
m_thpcRPC->nextDetectorElement(i, e)) {
512
513
514 std::map<Identifier, SimDataContent> channelSimDataMap;
515
516
517 while (i != e) {
519
520 TimedHitPtr<RPCSimHit> phit(*i++);
521
522
523 const RPCSimHit& hit(*phit);
524
525 const int idHit = hit.RPCid();
526
527 const double globalHitTime{
hitTime(phit)};
528
529 const double G4Time{hit.globalTime()};
530
531 const double bunchTime{globalHitTime - hit.globalTime()};
532
533 ATH_MSG_DEBUG(
"Global time " << globalHitTime <<
" G4 time " << G4Time <<
" Bunch time " << bunchTime);
534
536 ATH_MSG_VERBOSE(
"Validation: globalHitTime, G4Time, BCtime = " << globalHitTime <<
" " << G4Time <<
" " << bunchTime);
537 inputSimHitColl->Emplace(idHit, globalHitTime, hit.localPosition(),
539 hit.postLocalPosition(),
540 hit.energyDeposit(), hit.stepLength(), hit.particleEncoding(), hit.kineticEnergy());
541 }
542
543
551
553
554
556 const Identifier elementID =
m_idHelper->elementID(stationName,stationEta,stationPhi,doubletR,
isValid);
558 ATH_MSG_WARNING(
"Failed to construct the element ID from "<<stationName
559 <<", stationEta: "<<stationEta<<", stationPhi: "<<stationPhi<<", doubletR: "<<doubletR);
560 continue;
561 }
562
564 << " stationName " << stationName << " stationEta " << stationEta << " stationPhi " << stationPhi << " doubletR "
565 << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi << " gasGap " << gasGap);
566 const Identifier detElId{
m_idHelper->channelID(elementID, doubletZ, doubletPhi, 1,0, 1,
isValid)};
568 continue;
569 }
574 }
575
576
577 bool isValidEta{false}, isValidPhi{false};
578 const Identifier idpaneleta =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 0, 1, isValidEta);
579 const Identifier idpanelphi =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 1, 1, isValidPhi);
580 if (!isValidEta || !isValidPhi) {
582 << " stationName " << stationName << " stationEta " << stationEta << " stationPhi " << stationPhi
583 << " doubletR " << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi << " gasGap "
584 << gasGap);
585 continue;
586 }
587
588
592 const double corrtimejitter = tmp_CorrJitter > 0.01 ?
593 CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0., tmp_CorrJitter) : 0.;
594
595
596
597
598
599 const Amg::Vector3D hitDir{(hit.postLocalPosition() - hit.localPosition()).unit()};
601 Amg::intersect<3>(hit.localPosition(), hitDir, Amg::Vector3D::UnitX(), 0).value_or(0) * hitDir;
602
603 std::array<int, 3> pcseta =
physicalClusterSize(ctx, reEle, idpaneleta, gapCentre, rndmEngine);
604 ATH_MSG_VERBOSE(
"Simulated cluster on eta panel: size/first/last= " << pcseta[0] <<
"/" << pcseta[1] <<
"/" << pcseta[2]);
605 std::array<int, 3> pcsphi =
physicalClusterSize(ctx, reEle, idpanelphi, gapCentre, rndmEngine);
606 ATH_MSG_VERBOSE(
"Simulated cluster on phi panel: size/first/last= " << pcsphi[0] <<
"/" << pcsphi[1] <<
"/" << pcsphi[2]);
607
608
609
610
611 const Identifier atlasRpcIdeta =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 0, pcseta[1], isValidEta);
612 const Identifier atlasRpcIdphi =
m_idHelper->channelID(elementID, doubletZ, doubletPhi, gasGap, 1, pcsphi[1], isValidPhi);
613
615 const auto [etaStripOn, phiStripOn] =
detectionEfficiency(ctx, idpaneleta, idpanelphi, rndmEngine, particleLink);
616 ATH_MSG_DEBUG(
"SetPhiOn " << phiStripOn <<
" SetEtaOn " << etaStripOn);
617
618 for (bool imeasphi : {false, true}) {
619 if (!imeasphi && (!etaStripOn || !isValidEta)) continue;
620 if (imeasphi && (!phiStripOn || !isValidPhi)) continue;
621
622
623
624 const Identifier& atlasId = !imeasphi ? atlasRpcIdeta : atlasRpcIdphi;
625 std::array<int, 3>
pcs{!imeasphi ? pcseta : pcsphi};
626
627 ATH_MSG_DEBUG(
"SetOn: stationName " << stationName <<
" stationEta " << stationEta <<
" stationPhi " << stationPhi
628 << " doubletR " << doubletR << " doubletZ " << doubletZ << " doubletPhi " << doubletPhi
629 << " gasGap " << gasGap << " measphi " << imeasphi);
630
631
633 if (pcs[2] < 0){
634 continue;
635 }
636
637 ATH_MSG_DEBUG(
"Simulated cluster1: size/first/last= " << pcs[0] <<
"/" << pcs[1] <<
"/" << pcs[2]);
638
639
642
644 <<
" hit "<<
m_idHelper->print_to_string(atlasId)
652
653
655
656 double tns = G4Time + proptime + corrtimejitter;
657 ATH_MSG_VERBOSE(
"TOF+propagation time " << tns <<
" /s where proptime " << proptime <<
"/s");
658
659 double time = tns + bunchTime;
661
662
664
665 double*
b =
reinterpret_cast<double*
>(&packedMCword);
666
668
669
670
672
673
674
676 if (std::abs(hit.particleEncoding()) == 13 || hit.particleEncoding() == 0) {
677 if (channelSimDataMap.find(atlasId) == channelSimDataMap.end()) {
678 SimDataContent&
content = channelSimDataMap[atlasId];
680 content.deposits.push_back(deposit);
685 }
686 }
687 }
688
689
690
691
692
693
694
696 for (
int clus = pcs[1]; clus <=
pcs[2]; ++clus) {
697 Identifier newId =
m_idHelper->channelID(stationName, stationEta, stationPhi, doubletR, doubletZ,
698 doubletPhi, gasGap, imeasphi, clus,
isValid);
700 ATH_MSG_WARNING(__FILE__<<
":"<<__LINE__<<
"Channel "<< stationName<<
" "<<stationEta<<
" "<<stationPhi<<
" "<< doubletR<<
" "<<doubletZ
701 <<" "<< doubletPhi<<" "<< gasGap <<" "<< imeasphi<<" "<< clus<<" is invalid");
702 continue;
703 }
704
707 ATH_MSG_WARNING(
"Temporary skipping creation of RPC digit for stationName="
708 << stationName << ", eta=" << stationEta << ", phi=" << stationPhi << ", doubletR=" << doubletR
709 << ", doubletZ=" << doubletZ << ", doubletPhi=" << doubletPhi << ", gasGap=" << gasGap
710 << ", measuresPhi=" << imeasphi << ", strip=" << clus << ", cf. ATLASRECTS-6124");
711 return StatusCode::SUCCESS;
712 } else {
715 return StatusCode::FAILURE;
716 }
717 }
718
722
724 std::vector<MuonSimData::Deposit> newdeps;
725 newdeps.push_back(deposit);
726 m_sdo_tmp_map.insert(std::map<Identifier, std::vector<MuonSimData::Deposit>>::value_type(newId, newdeps));
727 } else {
729 }
730 }
731 }
732 }
733
735 for (auto it = channelSimDataMap.begin(); it != channelSimDataMap.end(); ++it) {
736 MuonSimData
simData(
it->second.deposits, 0);
739 auto insertResult = sdoContainer->insert(std::make_pair(
it->first,
simData));
740 if (!insertResult.second)
741 ATH_MSG_WARNING(
"Attention: this sdo is not recorded, since the identifier already exists in the sdoContainer map");
742 }
743 }
744
745 }
746
748
749 std::map<Identifier, std::vector<MuonSimData::Deposit>>::iterator map_iter =
m_sdo_tmp_map.begin();
751
753
754 const Identifier theId = (*map_iter).first;
756
757 const std::vector<MuonSimData::Deposit> theDeps = (*map_iter).second;
758
759
761 std::multimap<double, MuonSimData::Deposit>
times;
762
763
764 for (
unsigned int k = 0;
k < theDeps.size();
k++) {
765 double time = theDeps[
k].second.secondEntry();
766 times.insert(std::multimap<double, MuonSimData::Deposit>::value_type(time, theDeps[k]));
767 }
768
769
770
771 IdContext rpcContext =
m_idHelper->module_context();
772
773 std::multimap<double, MuonSimData::Deposit>::iterator map_dep_iter =
times.begin();
774
775
776 double last_time = -10000;
777 for (; map_dep_iter !=
times.end(); ++map_dep_iter) {
778 double currTime = (*map_dep_iter).first;
780
782
783
784 if (sdoContainer->find(theId) != sdoContainer->end())
785 {
786 std::map<Identifier, MuonSimData>::const_iterator
it = sdoContainer->find(theId);
787 std::vector<MuonSimData::Deposit> deps = ((*it).second).getdeposits();
788 deps.push_back((*map_dep_iter).second);
789 } else
790 {
791 std::vector<MuonSimData::Deposit> deposits;
792 deposits.push_back((*map_dep_iter).second);
793 std::pair<std::map<Identifier, MuonSimData>::iterator, bool> insertResult =
794 sdoContainer->insert(std::make_pair(theId, MuonSimData(deposits, 0)));
795 if (!insertResult.second)
797 "Attention TEMP: this sdo is not recorded, since the identifier already exists in the sdoContainer map");
798 }
799 }
800
801 if (std::abs(currTime - last_time) > (
m_deadTime)) {
802 ATH_MSG_DEBUG(
"deposit with time " << currTime <<
" is distant enough from previous (if any) hit on teh same strip");
803 last_time = (*map_dep_iter).first;
804
805
806 double uncorrjitter = 0;
809 if (tmp_UncorrJitter > 0.01) uncorrjitter = CLHEP::RandGaussZiggurat::shoot(rndmEngine, 0., tmp_UncorrJitter);
810
811
812
813
817
819 double newDigit_time = currTime + uncorrjitter +
m_rpc_time_shift - tp - propTimeFromStripCenter;
820
821 double digi_ToT = -1.;
823
824 ATH_MSG_VERBOSE(
"last_time=currTime " << last_time <<
" jitter " << uncorrjitter <<
" TOFcorrection " << tp <<
" shift "
826
827
828 bool outsideDigitizationWindow =
outsideWindow(newDigit_time);
829 if (outsideDigitizationWindow) {
830 ATH_MSG_VERBOSE(
"hit outside digitization window - do not produce digits");
834
835 continue;
836 }
837
838
839 last_time = (*map_dep_iter).first;
840
841 std::unique_ptr<RpcDigit> newDigit = std::make_unique<RpcDigit>(theId, newDigit_time, digi_ToT, false);
842
843 Identifier elemId =
m_idHelper->elementID(theId);
844 RpcDigitCollection* digitCollection = nullptr;
845
846 IdentifierHash coll_hash;
847 if (
m_idHelper->get_hash(elemId, coll_hash, &rpcContext)) {
848 ATH_MSG_ERROR(
"Unable to get RPC hash id from RPC Digit collection "
849 <<
"context begin_index = " << rpcContext.
begin_index()
850 <<
" context end_index = " << rpcContext.
end_index() <<
" the identifier is ");
852 }
853
854
856
857
858 if (coll_hash >= collections.size()) {
859 collections.resize (coll_hash+1);
860 }
861 digitCollection = collections[coll_hash].
get();
862 if (!digitCollection) {
863 collections[coll_hash] = std::make_unique<RpcDigitCollection>(elemId, coll_hash);
864 digitCollection = collections[coll_hash].
get();
865 }
866 digitCollection->
push_back(std::move(newDigit));
867
869
870 if (sdoContainer->find(theId) != sdoContainer->end()) {
871 std::map<Identifier, MuonSimData>::const_iterator
it = sdoContainer->find(theId);
872 std::vector<MuonSimData::Deposit> deps = ((*it).second).getdeposits();
873 deps.push_back((*map_dep_iter).second);
874 } else {
875 std::vector<MuonSimData::Deposit> deposits;
876 deposits.push_back((*map_dep_iter).second);
877 std::pair<std::map<Identifier, MuonSimData>::iterator, bool> insertResult =
878 sdoContainer->insert(std::make_pair(theId, MuonSimData(deposits, 0)));
879 if (!insertResult.second)
881 "Attention: this sdo is not recorded, since teh identifier already exists in the sdoContainer map");
882 }
883 }
884
885 } else
886 ATH_MSG_DEBUG(
"discarding digit due to dead time: " << (*map_dep_iter).first <<
" " << last_time);
887 }
888
889 }
890
891
896 }
897
898 return StatusCode::SUCCESS;
899}
float hitTime(const AFP_SIDSimHit &hit)
#define ATH_CHECK
Evaluate an expression and check for errors.
bool isValid(const T &p)
Av: we implement here an ATLAS-sepcific convention: all particles which are 99xxxxx are fine.
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
void show(std::ostream &out=std::cout) const
Print out in hex form.
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