364 edm::LogInfo(
"RPCRecHitValid") <<
"Cannot find simHit collection\n";
371 edm::LogInfo(
"RPCRecHitValid") <<
"Cannot find recHit collection\n";
378 edm::LogInfo(
"RPCRecHitValid") <<
"Cannot find TrackingParticle collection\n";
385 edm::LogInfo(
"RPCRecHitValid") <<
"Cannot find TrackingParticle to SimHit association map\n";
392 edm::LogInfo(
"RPCRecHitValid") <<
"Cannot find muon collection\n";
396 typedef edm::PSimHitContainer::const_iterator SimHitIter;
398 typedef std::vector<TrackPSimHitRef> SimHitRefs;
401 SimHitRefs muonSimHits, pthrSimHits;
403 for (
int i = 0,
n = simParticleHandle->size();
i <
n; ++
i) {
405 if (simParticle->pt() < 1.0
or simParticle->p() < 2.5)
409 SimHitRefs simHitsFromParticle;
410 auto range = std::equal_range(simHitsTPAssoc->begin(),
411 simHitsTPAssoc->end(),
414 for (
auto simParticleToHit =
range.first; simParticleToHit !=
range.second; ++simParticleToHit) {
415 auto simHit = simParticleToHit->second;
420 simHitsFromParticle.push_back(simParticleToHit->second);
422 const int nRPCHit = simHitsFromParticle.size();
423 const bool hasRPCHit = nRPCHit > 0;
425 if (
abs(simParticle->pdgId()) == 13) {
426 muonSimHits.insert(muonSimHits.end(), simHitsFromParticle.begin(), simHitsFromParticle.end());
429 int nRPCHitBarrel = 0;
430 int nRPCHitEndcap = 0;
431 for (
const auto &
simHit : simHitsFromParticle) {
432 const RPCDetId rpcDetId = static_cast<const RPCDetId>(
simHit->detUnitId());
433 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(rpcDetId));
437 if (rpcDetId.
region() == 0)
445 if (nRPCHitBarrel and nRPCHitEndcap) {
450 }
else if (nRPCHitBarrel) {
455 }
else if (nRPCHitEndcap) {
466 pthrSimHits.insert(pthrSimHits.end(), simHitsFromParticle.begin(), simHitsFromParticle.end());
470 switch (simParticle->pdgId()) {
509 int nRefHitBarrel = 0, nRefHitEndcap = 0;
510 for (
const auto &
simHit : muonSimHits) {
511 const RPCDetId detId = static_cast<const RPCDetId>(
simHit->detUnitId());
512 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(detId));
514 const int region = roll->id().region();
515 const int ring = roll->id().ring();
517 const int station = roll->id().station();
539 for (
const auto &
simHit : pthrSimHits) {
540 const RPCDetId detId = static_cast<const RPCDetId>(
simHit->detUnitId());
541 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(detId()));
543 const int region = roll->id().region();
544 const int ring = roll->id().ring();
546 const int station = roll->id().station();
569 int sumClusterSizeBarrel = 0, sumClusterSizeEndcap = 0;
570 int nRecHitBarrel = 0, nRecHitEndcap = 0;
571 for (RecHitIter recHitIter = recHitHandle->begin(); recHitIter != recHitHandle->end(); ++recHitIter) {
572 const RPCDetId detId = static_cast<const RPCDetId>(recHitIter->rpcId());
573 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(detId()));
577 const int region = roll->id().region();
578 const int ring = roll->id().ring();
580 const int station = roll->id().station();
584 const double time = recHitIter->timeError() >= 0 ? recHitIter->time() : recHitIter->BunchX() * 25;
590 sumClusterSizeBarrel += recHitIter->clusterSize();
599 sumClusterSizeEndcap += recHitIter->clusterSize();
607 if (roll->isIRPC()) {
613 const double nRecHit = nRecHitBarrel + nRecHitEndcap;
617 const int sumClusterSize = sumClusterSizeBarrel + sumClusterSizeEndcap;
620 if (nRecHitBarrel > 0) {
623 if (nRecHitEndcap > 0) {
629 typedef std::map<TrackPSimHitRef, RecHitIter> SimToRecHitMap;
630 SimToRecHitMap simToRecHitMap;
632 for (
const auto &
simHit : muonSimHits) {
633 const RPCDetId simDetId = static_cast<const RPCDetId>(
simHit->detUnitId());
637 const double simX =
simHit->localPosition().x();
639 for (RecHitIter recHitIter = recHitHandle->begin(); recHitIter != recHitHandle->end(); ++recHitIter) {
640 const RPCDetId recDetId = static_cast<const RPCDetId>(recHitIter->rpcId());
641 const RPCRoll *recRoll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(recDetId));
645 if (simDetId != recDetId)
648 const double recX = recHitIter->localPosition().x();
649 const double newDx = fabs(recX - simX);
652 SimToRecHitMap::const_iterator prevSimToReco = simToRecHitMap.find(
simHit);
653 if (prevSimToReco == simToRecHitMap.end()) {
654 simToRecHitMap.insert(std::make_pair(
simHit, recHitIter));
656 const double oldDx = fabs(prevSimToReco->second->localPosition().x() - simX);
659 simToRecHitMap[
simHit] = recHitIter;
667 int nMatchHitBarrel = 0, nMatchHitEndcap = 0;
668 for (
const auto &
match : simToRecHitMap) {
670 RecHitIter recHitIter =
match.second;
672 const RPCDetId detId = static_cast<const RPCDetId>(
simHit->detUnitId());
673 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(detId));
675 const int region = roll->id().region();
676 const int ring = roll->id().ring();
678 const int station = roll->id().station();
682 const double simX =
simHit->localPosition().x();
683 const double recX = recHitIter->localPosition().x();
684 const double errX =
sqrt(recHitIter->localPositionError().xx());
685 const double dX = recX - simX;
686 const double pull =
errX == 0 ? -999 : dX /
errX;
724 for (reco::MuonCollection::const_iterator
muon = muonHandle->begin();
muon != muonHandle->end(); ++
muon) {
725 if (!
muon->isGlobalMuon())
728 int nRPCHitBarrel = 0;
729 int nRPCHitEndcap = 0;
733 if (!(*recHit)->isValid())
735 const DetId detId = (*recHit)->geographicalId();
738 const RPCDetId rpcDetId = static_cast<const RPCDetId>(detId);
740 if (rpcDetId.
region() == 0)
746 const int nRPCHit = nRPCHitBarrel + nRPCHitEndcap;
748 if (nRPCHitBarrel and nRPCHitEndcap) {
753 }
else if (nRPCHitBarrel) {
758 }
else if (nRPCHitEndcap) {
771 for (RecHitIter recHitIter = recHitHandle->begin(); recHitIter != recHitHandle->end(); ++recHitIter) {
772 const RPCDetId detId = static_cast<const RPCDetId>(recHitIter->rpcId());
773 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(detId));
775 const int region = roll->id().region();
776 const int ring = roll->id().ring();
778 const int station = roll->id().station();
783 for (
const auto &
match : simToRecHitMap) {
784 if (recHitIter ==
match.second) {
791 int nPunchMatched = 0;
793 for (
const auto &
simHit : pthrSimHits) {
794 const int absSimHitPType =
abs(
simHit->particleType());
795 if (absSimHitPType == 13)
798 const RPCDetId simDetId = static_cast<const RPCDetId>(
simHit->detUnitId());
799 if (simDetId == detId)
803 if (nPunchMatched > 0) {
817 for (RecHitIter recHitIter = recHitHandle->begin(); recHitIter != recHitHandle->end(); ++recHitIter) {
818 const RPCDetId recDetId = static_cast<const RPCDetId>(recHitIter->rpcId());
819 const RPCRoll *roll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(recDetId));
821 const int region = roll->id().region();
828 const double recX = recHitIter->localPosition().x();
829 const double recErrX =
sqrt(recHitIter->localPositionError().xx());
832 for (SimHitIter simHitIter = simHitHandle->begin(); simHitIter != simHitHandle->end(); ++simHitIter) {
833 const RPCDetId simDetId = static_cast<const RPCDetId>(simHitIter->detUnitId());
834 const RPCRoll *simRoll = dynamic_cast<const RPCRoll *>(rpcGeom->
roll(simDetId));
838 if (simDetId != recDetId)
841 const double simX = simHitIter->localPosition().x();
842 const double dX = fabs(recX - simX);
844 if (dX / recErrX < 5) {