File indexing completed on 2023-01-27 23:42:31
0001
0002
0003 #include <memory>
0004 #include <string>
0005 #include <vector>
0006
0007 #include "SimTracker/TrackerHitAssociation/interface/TrackerHitAssociator.h"
0008
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011
0012 #include "DataFormats/Common/interface/Ref.h"
0013 #include "DataFormats/DetId/interface/DetId.h"
0014 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0015 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0016
0017 #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
0018
0019
0020 #include <numeric>
0021 #include <iostream>
0022
0023 using namespace std;
0024 using namespace edm;
0025
0026
0027
0028
0029 TrackerHitAssociator::Config::Config(edm::ConsumesCollector&& iC)
0030 : doPixel_(true), doStrip_(true), useOTph2_(false), doTrackAssoc_(false), assocHitbySimTrack_(false) {
0031 if (doStrip_) {
0032 if (useOTph2_)
0033 ph2OTrToken_ = iC.consumes<edm::DetSetVector<PixelDigiSimLink> >(edm::InputTag("simSiPixelDigis", "Tracker"));
0034 else
0035 stripToken_ = iC.consumes<edm::DetSetVector<StripDigiSimLink> >(edm::InputTag("simSiStripDigis"));
0036 }
0037 if (doPixel_) {
0038 if (useOTph2_)
0039 pixelToken_ = iC.consumes<edm::DetSetVector<PixelDigiSimLink> >(edm::InputTag("simSiPixelDigis", "Pixel"));
0040 else
0041 pixelToken_ = iC.consumes<edm::DetSetVector<PixelDigiSimLink> >(edm::InputTag("simSiPixelDigis"));
0042 }
0043 if (!doTrackAssoc_) {
0044 std::vector<std::string> trackerContainers;
0045 trackerContainers.reserve(12);
0046 trackerContainers.emplace_back("g4SimHitsTrackerHitsTIBLowTof");
0047 trackerContainers.emplace_back("g4SimHitsTrackerHitsTIBHighTof");
0048 trackerContainers.emplace_back("g4SimHitsTrackerHitsTIDLowTof");
0049 trackerContainers.emplace_back("g4SimHitsTrackerHitsTIDHighTof");
0050 trackerContainers.emplace_back("g4SimHitsTrackerHitsTOBLowTof");
0051 trackerContainers.emplace_back("g4SimHitsTrackerHitsTOBHighTof");
0052 trackerContainers.emplace_back("g4SimHitsTrackerHitsTECLowTof");
0053 trackerContainers.emplace_back("g4SimHitsTrackerHitsTECHighTof");
0054 trackerContainers.emplace_back("g4SimHitsTrackerHitsPixelBarrelLowTof");
0055 trackerContainers.emplace_back("g4SimHitsTrackerHitsPixelBarrelHighTof");
0056 trackerContainers.emplace_back("g4SimHitsTrackerHitsPixelEndcapLowTof");
0057 trackerContainers.emplace_back("g4SimHitsTrackerHitsPixelEndcapHighTof");
0058 cfTokens_.reserve(trackerContainers.size());
0059 simHitTokens_.reserve(trackerContainers.size());
0060 for (auto const& trackerContainer : trackerContainers) {
0061 cfTokens_.push_back(iC.consumes<CrossingFrame<PSimHit> >(edm::InputTag("mix", trackerContainer)));
0062 simHitTokens_.push_back(iC.consumes<std::vector<PSimHit> >(edm::InputTag("g4SimHits", trackerContainer)));
0063 }
0064 }
0065 }
0066
0067
0068
0069
0070 TrackerHitAssociator::Config::Config(const edm::ParameterSet& conf, edm::ConsumesCollector&& iC)
0071 : doPixel_(conf.getParameter<bool>("associatePixel")),
0072 doStrip_(conf.getParameter<bool>("associateStrip")),
0073 useOTph2_(conf.existsAs<bool>("usePhase2Tracker") ? conf.getParameter<bool>("usePhase2Tracker") : false),
0074
0075 doTrackAssoc_(conf.getParameter<bool>("associateRecoTracks")),
0076 assocHitbySimTrack_(
0077 conf.existsAs<bool>("associateHitbySimTrack") ? conf.getParameter<bool>("associateHitbySimTrack") : false) {
0078 if (doStrip_) {
0079 if (useOTph2_)
0080 ph2OTrToken_ = iC.consumes<edm::DetSetVector<PixelDigiSimLink> >(
0081 conf.getParameter<edm::InputTag>("phase2TrackerSimLinkSrc"));
0082 else
0083 stripToken_ =
0084 iC.consumes<edm::DetSetVector<StripDigiSimLink> >(conf.getParameter<edm::InputTag>("stripSimLinkSrc"));
0085 }
0086 if (doPixel_)
0087 pixelToken_ =
0088 iC.consumes<edm::DetSetVector<PixelDigiSimLink> >(conf.getParameter<edm::InputTag>("pixelSimLinkSrc"));
0089 if (!doTrackAssoc_) {
0090 std::vector<std::string> trackerContainers(conf.getParameter<std::vector<std::string> >("ROUList"));
0091 cfTokens_.reserve(trackerContainers.size());
0092 simHitTokens_.reserve(trackerContainers.size());
0093 for (auto const& trackerContainer : trackerContainers) {
0094 cfTokens_.push_back(iC.consumes<CrossingFrame<PSimHit> >(edm::InputTag("mix", trackerContainer)));
0095 simHitTokens_.push_back(iC.consumes<std::vector<PSimHit> >(edm::InputTag("g4SimHits", trackerContainer)));
0096 }
0097 }
0098 }
0099
0100
0101
0102
0103 TrackerHitAssociator::TrackerHitAssociator(const edm::Event& e, const TrackerHitAssociator::Config& config)
0104 : doPixel_(config.doPixel_),
0105 doStrip_(config.doStrip_),
0106 useOTph2_(config.useOTph2_),
0107 doTrackAssoc_(config.doTrackAssoc_),
0108 assocHitbySimTrack_(config.assocHitbySimTrack_) {
0109
0110 if (!doTrackAssoc_) {
0111 makeMaps(e, config);
0112 }
0113
0114 if (doStrip_) {
0115 if (useOTph2_)
0116 e.getByToken(config.ph2OTrToken_, ph2trackerdigisimlink);
0117 else
0118 e.getByToken(config.stripToken_, stripdigisimlink);
0119 }
0120 if (doPixel_)
0121 e.getByToken(config.pixelToken_, pixeldigisimlink);
0122 }
0123
0124 void TrackerHitAssociator::makeMaps(const edm::Event& theEvent, const TrackerHitAssociator::Config& config) {
0125
0126
0127
0128
0129 if (assocHitbySimTrack_) {
0130 for (auto const& cfToken : config.cfTokens_) {
0131 edm::Handle<CrossingFrame<PSimHit> > cf_simhit;
0132 int Nhits = 0;
0133 if (theEvent.getByToken(cfToken, cf_simhit)) {
0134 std::unique_ptr<MixCollection<PSimHit> > thisContainerHits(new MixCollection<PSimHit>(cf_simhit.product()));
0135 for (auto const& isim : *thisContainerHits) {
0136 DetId theDet(isim.detUnitId());
0137 SimHitMap[theDet].push_back(isim);
0138 ++Nhits;
0139 }
0140 LogDebug("TrkHitAssocTrace") << "simHits from crossing frames; map size = " << SimHitMap.size()
0141 << ", Hit count = " << Nhits << std::endl;
0142 }
0143 }
0144 for (auto const& simHitToken : config.simHitTokens_) {
0145 edm::Handle<std::vector<PSimHit> > simHits;
0146 int Nhits = 0;
0147 if (theEvent.getByToken(simHitToken, simHits)) {
0148 for (auto const& isim : *simHits) {
0149 DetId theDet(isim.detUnitId());
0150 SimHitMap[theDet].push_back(isim);
0151 ++Nhits;
0152 }
0153 LogDebug("TrkHitAssocTrace") << "simHits from prompt collections; map size = " << SimHitMap.size()
0154 << ", Hit count = " << Nhits << std::endl;
0155 }
0156 }
0157 } else {
0158 const char* const highTag = "HighTof";
0159 unsigned int tofBin;
0160 edm::EDConsumerBase::Labels labels;
0161 subDetTofBin theSubDetTofBin;
0162 unsigned int collectionIndex = 0;
0163 for (auto const& cfToken : config.cfTokens_) {
0164 collectionIndex++;
0165 edm::Handle<CrossingFrame<PSimHit> > cf_simhit;
0166 int Nhits = 0;
0167 if (theEvent.getByToken(cfToken, cf_simhit)) {
0168 std::unique_ptr<MixCollection<PSimHit> > thisContainerHits(new MixCollection<PSimHit>(cf_simhit.product()));
0169 theEvent.labelsForToken(cfToken, labels);
0170 if (std::strstr(labels.productInstance, highTag) != nullptr) {
0171 tofBin = StripDigiSimLink::HighTof;
0172 } else {
0173 tofBin = StripDigiSimLink::LowTof;
0174 }
0175 for (auto const& isim : *thisContainerHits) {
0176 DetId theDet(isim.detUnitId());
0177 theSubDetTofBin = std::make_pair(theDet.subdetId(), tofBin);
0178 SimHitCollMap[theSubDetTofBin] = collectionIndex;
0179 SimHitMap[SimHitCollMap[theSubDetTofBin]].push_back(isim);
0180 ++Nhits;
0181 }
0182 LogDebug("TrkHitAssocTrace") << "simHits from crossing frames " << collectionIndex << ": " << Nhits
0183 << std::endl;
0184 }
0185 }
0186 collectionIndex = 0;
0187 for (auto const& simHitToken : config.simHitTokens_) {
0188 collectionIndex++;
0189 edm::Handle<std::vector<PSimHit> > simHits;
0190 int Nhits = 0;
0191 if (theEvent.getByToken(simHitToken, simHits)) {
0192 theEvent.labelsForToken(simHitToken, labels);
0193 if (std::strstr(labels.productInstance, highTag) != nullptr) {
0194 tofBin = StripDigiSimLink::HighTof;
0195 } else {
0196 tofBin = StripDigiSimLink::LowTof;
0197 }
0198 for (auto const& isim : *simHits) {
0199 DetId theDet(isim.detUnitId());
0200 theSubDetTofBin = std::make_pair(theDet.subdetId(), tofBin);
0201 SimHitCollMap[theSubDetTofBin] = collectionIndex;
0202 SimHitMap[SimHitCollMap[theSubDetTofBin]].push_back(isim);
0203 ++Nhits;
0204 }
0205 LogDebug("TrkHitAssocTrace") << "simHits from prompt collection " << collectionIndex << ": " << Nhits
0206 << std::endl;
0207 }
0208 }
0209 }
0210 }
0211
0212 std::vector<PSimHit> TrackerHitAssociator::associateHit(const TrackingRecHit& thit) const {
0213 if (const SiTrackerMultiRecHit* rechit = dynamic_cast<const SiTrackerMultiRecHit*>(&thit)) {
0214 return associateMultiRecHit(rechit);
0215 }
0216
0217
0218 std::vector<PSimHit> result;
0219
0220 if (doTrackAssoc_)
0221 return result;
0222
0223
0224 std::vector<SimHitIdpr> simtrackid;
0225 std::vector<simhitAddr> simhitCFPos;
0226
0227
0228 DetId detid = thit.geographicalId();
0229 uint32_t detID = detid.rawId();
0230
0231
0232 associateHitId(thit, simtrackid, &simhitCFPos);
0233 LogDebug("TrkHitAssocTrace") << printDetBnchEvtTrk(detid, detID, simtrackid);
0234
0235
0236 if (!assocHitbySimTrack_ && !simhitCFPos.empty()) {
0237
0238
0239
0240
0241
0242
0243 if (dynamic_cast<const SiStripMatchedRecHit2D*>(&thit)) {
0244 for (auto const& theSimHitAddr : simhitCFPos) {
0245 simHitCollectionID theSimHitCollID = theSimHitAddr.first;
0246 auto it = SimHitMap.find(theSimHitCollID);
0247 if (it != SimHitMap.end()) {
0248 unsigned int theSimHitIndex = theSimHitAddr.second;
0249 if (theSimHitIndex < (it->second).size()) {
0250 const PSimHit& theSimHit = (it->second)[theSimHitIndex];
0251
0252 unsigned int simHitid = theSimHit.trackId();
0253 EncodedEventId simHiteid = theSimHit.eventId();
0254 for (auto const& id : simtrackid) {
0255 if (simHitid == id.first && simHiteid == id.second) {
0256 result.push_back(theSimHit);
0257 }
0258 }
0259 LogDebug("TrkHitAssocTrace") << "by CFpos, simHit detId = " << theSimHit.detUnitId() << " address = ("
0260 << theSimHitAddr.first << ", " << theSimHitIndex
0261 << "), process = " << theSimHit.processType() << " ("
0262 << theSimHit.eventId().bunchCrossing() << ", " << theSimHit.eventId().event()
0263 << ", " << theSimHit.trackId() << ")" << std::endl;
0264 }
0265 }
0266 }
0267 } else {
0268 for (auto const& theSimHitAddr : simhitCFPos) {
0269 simHitCollectionID theSimHitCollID = theSimHitAddr.first;
0270 auto it = SimHitMap.find(theSimHitCollID);
0271 if (it != SimHitMap.end()) {
0272 unsigned int theSimHitIndex = theSimHitAddr.second;
0273 if (theSimHitIndex < (it->second).size()) {
0274 result.push_back((it->second)[theSimHitIndex]);
0275 LogDebug("TrkHitAssocTrace") << "by CFpos, simHit detId = " << (it->second)[theSimHitIndex].detUnitId()
0276 << " address = (" << theSimHitCollID << ", " << theSimHitIndex
0277 << "), process = " << (it->second)[theSimHitIndex].processType() << " ("
0278 << (it->second)[theSimHitIndex].eventId().bunchCrossing() << ", "
0279 << (it->second)[theSimHitIndex].eventId().event() << ", "
0280 << (it->second)[theSimHitIndex].trackId() << ")" << std::endl;
0281 }
0282 }
0283 }
0284 }
0285 return result;
0286 }
0287
0288
0289 auto it = SimHitMap.find(detID);
0290 if (it != SimHitMap.end()) {
0291 for (auto const& ihit : it->second) {
0292 unsigned int simHitid = ihit.trackId();
0293 EncodedEventId simHiteid = ihit.eventId();
0294 for (auto id : simtrackid) {
0295 if (simHitid == id.first && simHiteid == id.second) {
0296 result.push_back(ihit);
0297 LogDebug("TrkHitAssocTrace") << "by TrackID, simHit detId = " << ihit.detUnitId()
0298 << ", process = " << ihit.processType() << " (" << ihit.eventId().bunchCrossing()
0299 << ", " << ihit.eventId().event() << ", " << ihit.trackId() << ")" << std::endl;
0300 break;
0301 }
0302 }
0303 }
0304
0305 } else {
0306
0307 auto itrphi = SimHitMap.find(detID + 2);
0308 auto itster = SimHitMap.find(detID + 1);
0309 if (itrphi != SimHitMap.end() && itster != SimHitMap.end()) {
0310 std::vector<PSimHit> simHitVector = itrphi->second;
0311 simHitVector.insert(simHitVector.end(), (itster->second).begin(), (itster->second).end());
0312 for (auto const& ihit : simHitVector) {
0313 unsigned int simHitid = ihit.trackId();
0314 EncodedEventId simHiteid = ihit.eventId();
0315 for (auto const& id : simtrackid) {
0316 if (simHitid == id.first && simHiteid == id.second) {
0317 result.push_back(ihit);
0318 LogDebug("TrkHitAssocTrace") << "by TrackID, simHit detId = " << ihit.detUnitId()
0319 << ", process = " << ihit.processType() << " ("
0320 << ihit.eventId().bunchCrossing() << ", " << ihit.eventId().event() << ", "
0321 << ihit.trackId() << ")" << std::endl;
0322 break;
0323 }
0324 }
0325 }
0326 }
0327 }
0328
0329 return result;
0330 }
0331
0332 std::vector<SimHitIdpr> TrackerHitAssociator::associateHitId(const TrackingRecHit& thit) const {
0333 std::vector<SimHitIdpr> simhitid;
0334 associateHitId(thit, simhitid);
0335 return simhitid;
0336 }
0337
0338 void TrackerHitAssociator::associateHitId(const TrackingRecHit& thit,
0339 std::vector<SimHitIdpr>& simtkid,
0340 std::vector<simhitAddr>* simhitCFPos) const {
0341 simtkid.clear();
0342
0343 if (const SiTrackerMultiRecHit* rechit = dynamic_cast<const SiTrackerMultiRecHit*>(&thit))
0344 simtkid = associateMultiRecHitId(rechit, simhitCFPos);
0345
0346
0347 if (const SiStripRecHit2D* rechit = dynamic_cast<const SiStripRecHit2D*>(&thit))
0348 associateSiStripRecHit(rechit, simtkid, simhitCFPos);
0349
0350
0351 else if (const SiStripRecHit1D* rechit = dynamic_cast<const SiStripRecHit1D*>(&thit))
0352 associateSiStripRecHit(rechit, simtkid, simhitCFPos);
0353
0354
0355 else if (const SiStripMatchedRecHit2D* rechit = dynamic_cast<const SiStripMatchedRecHit2D*>(&thit))
0356 simtkid = associateMatchedRecHit(rechit, simhitCFPos);
0357
0358
0359 else if (const ProjectedSiStripRecHit2D* rechit = dynamic_cast<const ProjectedSiStripRecHit2D*>(&thit))
0360 simtkid = associateProjectedRecHit(rechit, simhitCFPos);
0361
0362
0363 else if (const Phase2TrackerRecHit1D* rechit = dynamic_cast<const Phase2TrackerRecHit1D*>(&thit))
0364 associatePhase2TrackerRecHit(rechit, simtkid, simhitCFPos);
0365
0366
0367 else if (const SiPixelRecHit* rechit = dynamic_cast<const SiPixelRecHit*>(&thit))
0368 associatePixelRecHit(rechit, simtkid, simhitCFPos);
0369
0370
0371 if (trackerHitRTTI::isFast(thit))
0372 simtkid = associateFastRecHit(static_cast<const FastTrackerRecHit*>(&thit));
0373 }
0374
0375 template <typename T>
0376 inline void TrackerHitAssociator::associateSiStripRecHit(const T* simplerechit,
0377 std::vector<SimHitIdpr>& simtrackid,
0378 std::vector<simhitAddr>* simhitCFPos) const {
0379 const SiStripCluster* clust = &(*simplerechit->cluster());
0380 associateSimpleRecHitCluster(clust, simplerechit->geographicalId(), simtrackid, simhitCFPos);
0381 }
0382
0383
0384
0385
0386 void TrackerHitAssociator::associateCluster(const SiStripCluster* clust,
0387 const DetId& detid,
0388 std::vector<SimHitIdpr>& simtrackid,
0389 std::vector<PSimHit>& simhit) const {
0390 std::vector<simhitAddr> simhitCFPos;
0391 associateSimpleRecHitCluster(clust, detid, simtrackid, &simhitCFPos);
0392
0393 for (auto const& theSimHitAddr : simhitCFPos) {
0394 simHitCollectionID theSimHitCollID = theSimHitAddr.first;
0395 auto it = SimHitMap.find(theSimHitCollID);
0396
0397 if (it != SimHitMap.end()) {
0398 unsigned int theSimHitIndex = theSimHitAddr.second;
0399 if (theSimHitIndex < (it->second).size())
0400 simhit.push_back((it->second)[theSimHitIndex]);
0401 LogDebug("TrkHitAssocTrace") << "For cluster, simHit detId = " << (it->second)[theSimHitIndex].detUnitId()
0402 << " address = (" << theSimHitCollID << ", " << theSimHitIndex
0403 << "), process = " << (it->second)[theSimHitIndex].processType()
0404 << " (bnch, evt, trk) = (" << (it->second)[theSimHitIndex].eventId().bunchCrossing()
0405 << ", " << (it->second)[theSimHitIndex].eventId().event() << ", "
0406 << (it->second)[theSimHitIndex].trackId() << ")" << std::endl;
0407 }
0408 }
0409 }
0410
0411 void TrackerHitAssociator::associateSimpleRecHitCluster(const SiStripCluster* clust,
0412 const DetId& detid,
0413 std::vector<SimHitIdpr>& simtrackid,
0414 std::vector<simhitAddr>* simhitCFPos) const {
0415 uint32_t detID = detid.rawId();
0416 auto isearch = stripdigisimlink->find(detID);
0417 if (isearch != stripdigisimlink->end()) {
0418 auto link_detset = (*isearch);
0419
0420 if (clust != nullptr) {
0421 int clusiz = clust->amplitudes().size();
0422 int first = clust->firstStrip();
0423 int last = first + clusiz;
0424
0425 LogDebug("TrkHitAssocDbg") << "Cluster size " << clusiz << " first strip = " << first
0426 << " last strip = " << last - 1 << std::endl
0427 << " detID = " << detID << " DETSET size = " << link_detset.data.size() << std::endl;
0428 int channel;
0429 for (const auto& linkiter : link_detset.data) {
0430 channel = (int)(linkiter.channel());
0431 if (channel >= first && channel < last) {
0432 LogDebug("TrkHitAssocDbg") << "Channel = " << std::setw(4) << linkiter.channel()
0433 << ", TrackID = " << std::setw(8) << linkiter.SimTrackId()
0434 << ", tofBin = " << std::setw(3) << linkiter.TofBin()
0435 << ", fraction = " << std::setw(8) << linkiter.fraction()
0436 << ", Position = " << linkiter.CFposition() << std::endl;
0437 SimHitIdpr currentId(linkiter.SimTrackId(), linkiter.eventId());
0438
0439
0440
0441
0442 if (find(simtrackid.begin(), simtrackid.end(), currentId) == simtrackid.end()) {
0443 LogDebug("TrkHitAssocDbg") << " Adding track id = " << currentId.first
0444 << " Event id = " << currentId.second.event()
0445 << " Bunch Xing = " << currentId.second.bunchCrossing() << std::endl;
0446 simtrackid.push_back(currentId);
0447 }
0448
0449 if (simhitCFPos != nullptr) {
0450
0451
0452 unsigned int currentCFPos = linkiter.CFposition();
0453 unsigned int tofBin = linkiter.TofBin();
0454 subDetTofBin theSubDetTofBin = std::make_pair(detid.subdetId(), tofBin);
0455 auto it = SimHitCollMap.find(theSubDetTofBin);
0456 if (it != SimHitCollMap.end()) {
0457 simhitAddr currentAddr = std::make_pair(it->second, currentCFPos);
0458 if (find(simhitCFPos->begin(), simhitCFPos->end(), currentAddr) == simhitCFPos->end()) {
0459 simhitCFPos->push_back(currentAddr);
0460 }
0461 }
0462 }
0463 }
0464 }
0465 } else {
0466 edm::LogError("TrackerHitAssociator") << "no cluster reference attached";
0467 }
0468 }
0469 }
0470
0471 std::vector<SimHitIdpr> TrackerHitAssociator::associateMatchedRecHit(const SiStripMatchedRecHit2D* matchedrechit,
0472 std::vector<simhitAddr>* simhitCFPos) const {
0473 std::vector<SimHitIdpr> matched_mono;
0474 std::vector<SimHitIdpr> matched_st;
0475
0476 const SiStripRecHit2D mono = matchedrechit->monoHit();
0477 const SiStripRecHit2D st = matchedrechit->stereoHit();
0478
0479 associateSiStripRecHit(&mono, matched_mono, simhitCFPos);
0480 associateSiStripRecHit(&st, matched_st, simhitCFPos);
0481
0482
0483 std::vector<SimHitIdpr> simtrackid;
0484 if (!(matched_mono.empty() || matched_st.empty())) {
0485 for (auto const& mhit : matched_mono) {
0486
0487 if (find(simtrackid.begin(), simtrackid.end(), mhit) == simtrackid.end()) {
0488
0489 if (find(matched_st.begin(), matched_st.end(), mhit) != matched_st.end()) {
0490 simtrackid.push_back(mhit);
0491 }
0492 }
0493 }
0494 }
0495 return simtrackid;
0496 }
0497
0498 std::vector<SimHitIdpr> TrackerHitAssociator::associateProjectedRecHit(const ProjectedSiStripRecHit2D* projectedrechit,
0499 std::vector<simhitAddr>* simhitCFPos) const {
0500
0501
0502 std::vector<SimHitIdpr> matched_mono;
0503
0504 const SiStripRecHit2D mono = projectedrechit->originalHit();
0505 associateSiStripRecHit(&mono, matched_mono, simhitCFPos);
0506 return matched_mono;
0507 }
0508
0509 void TrackerHitAssociator::associatePhase2TrackerRecHit(const Phase2TrackerRecHit1D* rechit,
0510 std::vector<SimHitIdpr>& simtrackid,
0511 std::vector<simhitAddr>* simhitCFPos) const {
0512
0513
0514
0515 DetId detid = rechit->geographicalId();
0516 uint32_t detID = detid.rawId();
0517
0518 auto isearch = ph2trackerdigisimlink->find(detID);
0519 if (isearch != ph2trackerdigisimlink->end()) {
0520 auto link_detset = (*isearch);
0521 Phase2TrackerRecHit1D::ClusterRef const& cluster = rechit->cluster();
0522
0523
0524 if (!(cluster.isNull())) {
0525 int minRow = (*cluster).firstStrip();
0526 int maxRow = (*cluster).firstStrip() + (*cluster).size();
0527 int Col = (*cluster).column();
0528 LogDebug("TrkHitAssocDbg") << " Cluster minRow " << minRow << " maxRow " << maxRow << " column " << Col
0529 << std::endl;
0530 int dsl = 0;
0531 for (auto const& linkiter : link_detset.data) {
0532 ++dsl;
0533 std::pair<int, int> coord = Phase2TrackerDigi::channelToPixel(linkiter.channel());
0534 LogDebug("TrkHitAssocDbg") << " " << dsl << ") Digi link: row " << coord.first << " col " << coord.second
0535 << std::endl;
0536 if (coord.first <= maxRow && coord.first >= minRow && coord.second == Col) {
0537 LogDebug("TrkHitAssocDbg") << " !-> trackid " << linkiter.SimTrackId() << endl
0538 << " fraction " << linkiter.fraction() << endl;
0539 SimHitIdpr currentId(linkiter.SimTrackId(), linkiter.eventId());
0540 if (find(simtrackid.begin(), simtrackid.end(), currentId) == simtrackid.end()) {
0541 simtrackid.push_back(currentId);
0542 }
0543
0544 if (simhitCFPos != nullptr) {
0545
0546
0547 unsigned int currentCFPos = linkiter.CFposition();
0548 unsigned int tofBin = linkiter.TofBin();
0549 subDetTofBin theSubDetTofBin = std::make_pair(detid.subdetId(), tofBin);
0550 auto it = SimHitCollMap.find(theSubDetTofBin);
0551 if (it != SimHitCollMap.end()) {
0552 simhitAddr currentAddr = std::make_pair(it->second, currentCFPos);
0553 if (find(simhitCFPos->begin(), simhitCFPos->end(), currentAddr) == simhitCFPos->end()) {
0554 simhitCFPos->push_back(currentAddr);
0555 }
0556 }
0557 }
0558 }
0559 }
0560 } else {
0561 edm::LogError("TrackerHitAssociator") << "no Phase2 outer tracker cluster reference attached";
0562 }
0563 }
0564 }
0565
0566 void TrackerHitAssociator::associatePixelRecHit(const SiPixelRecHit* pixelrechit,
0567 std::vector<SimHitIdpr>& simtrackid,
0568 std::vector<simhitAddr>* simhitCFPos) const {
0569
0570
0571
0572 DetId detid = pixelrechit->geographicalId();
0573 uint32_t detID = detid.rawId();
0574
0575 auto isearch = pixeldigisimlink->find(detID);
0576 if (isearch != pixeldigisimlink->end()) {
0577 auto link_detset = (*isearch);
0578 SiPixelRecHit::ClusterRef const& cluster = pixelrechit->cluster();
0579
0580
0581
0582 if (!(cluster.isNull())) {
0583
0584 int minPixelRow = (*cluster).minPixelRow();
0585 int maxPixelRow = (*cluster).maxPixelRow();
0586 int minPixelCol = (*cluster).minPixelCol();
0587 int maxPixelCol = (*cluster).maxPixelCol();
0588 LogDebug("TrkHitAssocDbg") << " Cluster minRow " << minPixelRow << " maxRow " << maxPixelRow << std::endl
0589 << " Cluster minCol " << minPixelCol << " maxCol " << maxPixelCol << std::endl;
0590 int dsl = 0;
0591 for (auto const& linkiter : link_detset.data) {
0592 ++dsl;
0593 std::pair<int, int> pixel_coord = PixelDigi::channelToPixel(linkiter.channel());
0594 LogDebug("TrkHitAssocDbg") << " " << dsl << ") Digi link: row " << pixel_coord.first << " col "
0595 << pixel_coord.second << std::endl;
0596 if (pixel_coord.first <= maxPixelRow && pixel_coord.first >= minPixelRow && pixel_coord.second <= maxPixelCol &&
0597 pixel_coord.second >= minPixelCol) {
0598 LogDebug("TrkHitAssocDbg") << " !-> trackid " << linkiter.SimTrackId() << endl
0599 << " fraction " << linkiter.fraction() << endl;
0600 SimHitIdpr currentId(linkiter.SimTrackId(), linkiter.eventId());
0601 if (find(simtrackid.begin(), simtrackid.end(), currentId) == simtrackid.end()) {
0602 simtrackid.push_back(currentId);
0603 }
0604
0605 if (simhitCFPos != nullptr) {
0606
0607
0608 unsigned int currentCFPos = linkiter.CFposition();
0609 unsigned int tofBin = linkiter.TofBin();
0610 subDetTofBin theSubDetTofBin = std::make_pair(detid.subdetId(), tofBin);
0611 auto it = SimHitCollMap.find(theSubDetTofBin);
0612 if (it != SimHitCollMap.end()) {
0613 simhitAddr currentAddr = std::make_pair(it->second, currentCFPos);
0614 if (find(simhitCFPos->begin(), simhitCFPos->end(), currentAddr) == simhitCFPos->end()) {
0615 simhitCFPos->push_back(currentAddr);
0616 }
0617 }
0618 }
0619 }
0620 }
0621 } else {
0622 edm::LogError("TrackerHitAssociator") << "no Pixel cluster reference attached";
0623 }
0624 }
0625 }
0626
0627 std::vector<PSimHit> TrackerHitAssociator::associateMultiRecHit(const SiTrackerMultiRecHit* multirechit) const {
0628 std::vector<const TrackingRecHit*> componenthits = multirechit->recHits();
0629
0630 int size = multirechit->weights().size(), idmostprobable = 0;
0631
0632 for (int i = 0; i < size; ++i) {
0633 if (multirechit->weight(i) > multirechit->weight(idmostprobable))
0634 idmostprobable = i;
0635 }
0636
0637 return associateHit(*componenthits[idmostprobable]);
0638 }
0639
0640 std::vector<SimHitIdpr> TrackerHitAssociator::associateMultiRecHitId(const SiTrackerMultiRecHit* multirechit,
0641 std::vector<simhitAddr>* simhitCFPos) const {
0642 std::vector<const TrackingRecHit*> componenthits = multirechit->recHits();
0643 int size = multirechit->weights().size(), idmostprobable = 0;
0644
0645 for (int i = 0; i < size; ++i) {
0646 if (multirechit->weight(i) > multirechit->weight(idmostprobable))
0647 idmostprobable = i;
0648 }
0649
0650 std::vector<SimHitIdpr> simhitid;
0651 associateHitId(*componenthits[idmostprobable], simhitid, simhitCFPos);
0652 return simhitid;
0653 }
0654
0655
0656 std::vector<SimHitIdpr> TrackerHitAssociator::associateFastRecHit(const FastTrackerRecHit* rechit) const {
0657 vector<SimHitIdpr> simtrackid;
0658 simtrackid.clear();
0659 for (size_t index = 0, indexEnd = rechit->nSimTrackIds(); index < indexEnd; ++index) {
0660 SimHitIdpr currentId(rechit->simTrackId(index), EncodedEventId(rechit->simTrackEventId(index)));
0661 simtrackid.push_back(currentId);
0662 }
0663 return simtrackid;
0664 }
0665
0666 inline std::string TrackerHitAssociator::printDetBnchEvtTrk(const DetId& detid,
0667 const uint32_t& detID,
0668 std::vector<SimHitIdpr>& simtrackid) const {
0669 std::stringstream message;
0670 message << "recHit subdet, detID = " << detid.subdetId() << ", " << detID << ", (bnch, evt, trk) = ";
0671 for (size_t i = 0; i < simtrackid.size(); ++i)
0672 message << ", (" << simtrackid[i].second.bunchCrossing() << ", " << simtrackid[i].second.event() << ", "
0673 << simtrackid[i].first << ")";
0674
0675 return message.str();
0676 }