408 LogDebug(
"SiStripHitEfficiencyWorker") <<
"beginning analyze from HitEff";
424 }
else if (metaData.
isValid()) {
428 edm::LogWarning(
"SiStripHitEfficiencyWorker") <<
"could not find a source for the Luminosity and PU";
460 if (not fedErrorIdsCol_h.isValid() and not fedErrorIdsVec_h.isValid()) {
462 <<
"no valid product for SiStrip DetIds with FED errors (see parameter \"siStripDigis\"), " 463 "neither for new format (DetIdVector) nor old format (DetIdCollection)";
465 auto const& fedErrorIds = fedErrorIdsVec_h.isValid() ? *fedErrorIdsVec_h : fedErrorIdsCol_h->as_vector();
468 for (
const auto& fedErr : fedErrorIds) {
492 LogDebug(
"SiStripHitEfficiencyWorker") <<
"number ckf tracks found = " << tracksCKF->size();
502 if (!tracksCKF->empty()) {
506 LogDebug(
"SiStripHitEfficiencyWorker")
512 for (
const auto& trajTrack : *trajTrackAssociationHandle) {
516 auto TMeas = trajTrack.key->measurements();
526 bool hasMissingHits{
false};
527 int previous_layer{999};
528 std::vector<unsigned int> missedLayers;
530 for (
const auto& itm : TMeas) {
531 auto theHit = itm.recHit();
532 unsigned int iidd = theHit->geographicalId().rawId();
533 int layer = ::checkLayer(iidd, tTopo);
534 int missedLayer =
layer + 1;
535 int previousMissedLayer = (
layer + 2);
536 int diffPreviousLayer = (
layer - previous_layer);
539 if (diffPreviousLayer == -2 && missedLayer > k_LayersStart && missedLayer < k_LayersAtTOBEnd) {
541 hasMissingHits =
true;
544 else if (diffPreviousLayer == -2 && (missedLayer > k_LayersAtTOBEnd + 1 && missedLayer <= k_LayersAtTIDEnd)) {
546 hasMissingHits =
true;
549 else if (diffPreviousLayer == -2 && missedLayer > k_LayersAtTIDEnd && missedLayer <= k_LayersAtTECEnd) {
551 hasMissingHits =
true;
555 if ((
layer > k_LayersStart &&
layer <= k_LayersAtTIBEnd) && (previous_layer == 12)) {
557 hasMissingHits =
true;
561 if ((
layer > k_LayersAtTIBEnd &&
layer <= k_LayersAtTOBEnd) && (previous_layer == 15)) {
563 hasMissingHits =
true;
569 if (diffPreviousLayer == -3 && missedLayer > k_LayersStart && missedLayer < k_LayersAtTOBEnd &&
570 previousMissedLayer > k_LayersStart && previousMissedLayer < k_LayersAtTOBEnd) {
573 hasMissingHits =
true;
577 else if (diffPreviousLayer == -3 && missedLayer > k_LayersAtTIDEnd && missedLayer <= k_LayersAtTECEnd &&
578 previousMissedLayer > k_LayersAtTIDEnd && previousMissedLayer <= k_LayersAtTECEnd) {
581 hasMissingHits =
true;
585 hasMissingHits =
true;
588 missedLayers.push_back(
layer);
589 previous_layer =
layer;
594 unsigned int prev_TKlayers = 0;
595 for (
auto itm = TMeas.cbegin(); itm != TMeas.cend(); ++itm) {
596 const auto theInHit = (*itm).recHit();
601 LogDebug(
"SiStripHitEfficiencyWorker") <<
"theInHit is valid = " << theInHit->isValid();
603 unsigned int iidd = theInHit->geographicalId().rawId();
605 unsigned int TKlayers = ::checkLayer(iidd, tTopo);
607 bool foundConsMissingHits{
false};
613 LogDebug(
"SiStripHitEfficiencyWorker") <<
"TKlayer from trajectory: " << TKlayers <<
" from module = " << iidd
614 <<
" matched/stereo/rphi = " << ((iidd & 0x3) == 0) <<
"/" 615 << ((iidd & 0x3) == 1) <<
"/" << ((iidd & 0x3) == 2);
619 bool isFirstMeas = (itm == (TMeas.end() - 1));
620 bool isLastMeas = (itm == (TMeas.begin()));
633 if (TKlayers == bounds::k_LayersAtTOBEnd || TKlayers == bounds::k_LayersAtTECEnd) {
634 LogDebug(
"SiStripHitEfficiencyWorker") <<
"skipping original TM for TOB 6 or TEC 9";
638 std::vector<TrajectoryAtInvalidHit> TMs;
646 if (::isDoubleSided(iidd, tTopo) && ((iidd & 0x3) == 0)) {
649 TMs.emplace_back(*itm, tTopo, tkgeom,
propagator, 1);
650 TMs.emplace_back(*itm, tTopo, tkgeom,
propagator, 2);
651 }
else if (::isDoubleSided(iidd, tTopo) && (!::check2DPartner(iidd, TMeas))) {
654 TMs.emplace_back(*itm, tTopo, tkgeom,
propagator, 1);
655 TMs.emplace_back(*itm, tTopo, tkgeom,
propagator, 2);
656 LogDebug(
"SiStripHitEfficiencyWorker") <<
" found a hit with a missing partner";
659 TMs.emplace_back(*itm, tTopo, tkgeom,
propagator);
662 bool missingHitAdded{
false};
663 std::vector<TrajectoryMeasurement> tmpTmeas, prev_tmpTmeas;
664 unsigned int misLayer = TKlayers + 1;
665 unsigned int previousMisLayer = TKlayers + 2;
668 if (
int(TKlayers) -
int(prev_TKlayers) == -2) {
669 const DetLayer* detlayer = itm->layer();
672 std::vector<DetLayer::DetWithState> compatDets = detlayer->
compatibleDets(tsos, prop, chi2Estimator);
674 if (misLayer > k_LayersAtTIDEnd && misLayer < k_LayersAtTECEnd) {
675 std::vector<ForwardDetLayer const*> negTECLayers = measTracker.geometricSearchTracker()->negTecLayers();
676 std::vector<ForwardDetLayer const*> posTECLayers = measTracker.geometricSearchTracker()->posTecLayers();
677 const DetLayer* tecLayerneg = negTECLayers[misLayer - k_LayersAtTIDEnd - 1];
678 const DetLayer* tecLayerpos = posTECLayers[misLayer - k_LayersAtTIDEnd - 1];
679 if (tTopo->tecSide(iidd) == 1) {
680 tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, prop, chi2Estimator);
681 }
else if (tTopo->tecSide(iidd) == 2) {
682 tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, prop, chi2Estimator);
686 else if (misLayer == (k_LayersAtTIDEnd - 1) ||
687 misLayer == k_LayersAtTIDEnd) {
688 std::vector<ForwardDetLayer const*> negTIDLayers = measTracker.geometricSearchTracker()->negTidLayers();
689 std::vector<ForwardDetLayer const*> posTIDLayers = measTracker.geometricSearchTracker()->posTidLayers();
690 const DetLayer* tidLayerneg = negTIDLayers[misLayer - k_LayersAtTOBEnd - 1];
691 const DetLayer* tidLayerpos = posTIDLayers[misLayer - k_LayersAtTOBEnd - 1];
693 if (tTopo->tidSide(iidd) == 1) {
694 tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, prop, chi2Estimator);
695 }
else if (tTopo->tidSide(iidd) == 2) {
696 tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, prop, chi2Estimator);
700 if (misLayer > k_LayersStart && misLayer < k_LayersAtTOBEnd) {
702 std::vector<BarrelDetLayer const*> barrelTIBLayers = measTracker.geometricSearchTracker()->tibLayers();
703 std::vector<BarrelDetLayer const*> barrelTOBLayers = measTracker.geometricSearchTracker()->tobLayers();
705 if (misLayer > k_LayersStart && misLayer <= k_LayersAtTIBEnd) {
706 const DetLayer* tibLayer = barrelTIBLayers[misLayer - k_LayersStart - 1];
707 tmpTmeas = layerMeasurements.measurements(*tibLayer, tsos, prop, chi2Estimator);
708 }
else if (misLayer > k_LayersAtTIBEnd && misLayer < k_LayersAtTOBEnd) {
709 const DetLayer* tobLayer = barrelTOBLayers[misLayer - k_LayersAtTIBEnd - 1];
710 tmpTmeas = layerMeasurements.measurements(*tobLayer, tsos, prop, chi2Estimator);
714 if ((
int(TKlayers) > k_LayersStart &&
int(TKlayers) <= k_LayersAtTIBEnd) &&
int(prev_TKlayers) == 12) {
715 const DetLayer* detlayer = itm->layer();
718 std::vector<DetLayer::DetWithState> compatDets = detlayer->
compatibleDets(tsos, prop, chi2Estimator);
719 std::vector<ForwardDetLayer const*> negTIDLayers = measTracker.geometricSearchTracker()->negTidLayers();
720 std::vector<ForwardDetLayer const*> posTIDLayers = measTracker.geometricSearchTracker()->posTidLayers();
722 const DetLayer* tidLayerneg = negTIDLayers[k_LayersStart];
723 const DetLayer* tidLayerpos = posTIDLayers[k_LayersStart];
724 if (tTopo->tidSide(iidd) == 1) {
725 tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, prop, chi2Estimator);
726 }
else if (tTopo->tidSide(iidd) == 2) {
727 tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, prop, chi2Estimator);
731 if ((
int(TKlayers) > k_LayersAtTIBEnd &&
int(TKlayers) <= k_LayersAtTOBEnd) &&
int(prev_TKlayers) == 15) {
732 const DetLayer* detlayer = itm->layer();
735 std::vector<DetLayer::DetWithState> compatDets = detlayer->
compatibleDets(tsos, prop, chi2Estimator);
737 std::vector<ForwardDetLayer const*> negTECLayers = measTracker.geometricSearchTracker()->negTecLayers();
738 std::vector<ForwardDetLayer const*> posTECLayers = measTracker.geometricSearchTracker()->posTecLayers();
740 const DetLayer* tecLayerneg = negTECLayers[k_LayersStart];
741 const DetLayer* tecLayerpos = posTECLayers[k_LayersStart];
742 if (tTopo->tecSide(iidd) == 1) {
743 tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, prop, chi2Estimator);
744 }
else if (tTopo->tecSide(iidd) == 2) {
745 tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, prop, chi2Estimator);
750 if (
int(TKlayers) -
int(prev_TKlayers) == -3) {
751 foundConsMissingHits =
true;
752 const DetLayer* detlayer = itm->layer();
755 std::vector<DetLayer::DetWithState> compatDets = detlayer->
compatibleDets(tsos, prop, chi2Estimator);
757 if (misLayer > k_LayersStart && misLayer <= k_LayersAtTOBEnd && previousMisLayer > k_LayersStart &&
758 previousMisLayer <= k_LayersAtTOBEnd) {
759 std::vector<BarrelDetLayer const*> barrelTIBLayers = measTracker.geometricSearchTracker()->tibLayers();
760 std::vector<BarrelDetLayer const*> barrelTOBLayers = measTracker.geometricSearchTracker()->tobLayers();
761 if (misLayer > k_LayersStart && misLayer < k_LayersAtTIBEnd) {
762 const DetLayer* tibLayer = barrelTIBLayers[misLayer - k_LayersStart - 1];
763 const DetLayer* prevTibLayer = barrelTIBLayers[previousMisLayer - k_LayersStart - 1];
765 tmpTmeas = layerMeasurements.measurements(*tibLayer, tsos, prop, chi2Estimator);
766 prev_tmpTmeas = layerMeasurements.measurements(*prevTibLayer, tsos, prop, chi2Estimator);
767 }
else if (misLayer > k_LayersAtTIBEnd && misLayer < k_LayersAtTOBEnd) {
768 const DetLayer* tobLayer = barrelTOBLayers[misLayer - k_LayersAtTIBEnd - 1];
769 const DetLayer* prevTobLayer = barrelTOBLayers[previousMisLayer - k_LayersAtTIBEnd - 1];
770 tmpTmeas = layerMeasurements.measurements(*tobLayer, tsos, prop, chi2Estimator);
771 prev_tmpTmeas = layerMeasurements.measurements(*prevTobLayer, tsos, prop, chi2Estimator);
773 }
else if (misLayer > k_LayersAtTIDEnd && misLayer < k_LayersAtTECEnd &&
774 previousMisLayer > k_LayersAtTIDEnd && previousMisLayer < k_LayersAtTECEnd) {
775 std::vector<ForwardDetLayer const*> negTECLayers = measTracker.geometricSearchTracker()->negTecLayers();
776 std::vector<ForwardDetLayer const*> posTECLayers = measTracker.geometricSearchTracker()->posTecLayers();
778 const DetLayer* tecLayerneg = negTECLayers[misLayer - k_LayersAtTIDEnd - 1];
779 const DetLayer* prevTecLayerneg = negTECLayers[previousMisLayer - k_LayersAtTIDEnd - 1];
781 const DetLayer* tecLayerpos = posTECLayers[misLayer - k_LayersAtTIDEnd - 1];
782 const DetLayer* prevTecLayerpos = posTECLayers[previousMisLayer - k_LayersAtTIDEnd - 1];
784 if (tTopo->tecSide(iidd) == 1) {
785 tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, prop, chi2Estimator);
786 prev_tmpTmeas = layerMeasurements.measurements(*prevTecLayerneg, tsos, prop, chi2Estimator);
787 }
else if (tTopo->tecSide(iidd) == 2) {
788 tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, prop, chi2Estimator);
789 prev_tmpTmeas = layerMeasurements.measurements(*prevTecLayerpos, tsos, prop, chi2Estimator);
793 if (!tmpTmeas.empty() && !foundConsMissingHits) {
795 unsigned int iidd_tmp = TM_tmp.
recHit()->geographicalId().rawId();
797 LogDebug(
"SiStripHitEfficiency:HitEff") <<
" hit actually being added to TM vector";
800 if (::isDoubleSided(iidd_tmp, tTopo)) {
805 missingHitAdded =
true;
810 if (!tmpTmeas.empty() && !prev_tmpTmeas.empty() &&
811 foundConsMissingHits) {
815 unsigned int modIdInner = TM_tmp1.
recHit()->geographicalId().rawId();
816 unsigned int modIdOuter = TM_tmp2.recHit()->geographicalId().rawId();
817 bool innerModInactive =
false, outerModInactive =
false;
818 for (
const auto& tm : tmpTmeas) {
819 unsigned int tmModId = tm.recHit()->geographicalId().rawId();
820 if (tmModId == modIdInner && tm.recHit()->getType() == 2) {
821 innerModInactive =
true;
825 for (
const auto& tm : prev_tmpTmeas) {
826 unsigned int tmModId = tm.recHit()->geographicalId().rawId();
827 if (tmModId == modIdOuter && tm.recHit()->getType() == 2) {
828 outerModInactive =
true;
833 if (outerModInactive) {
834 if (modIdInner != 0) {
835 LogDebug(
"SiStripHitEfficiency:HitEff") <<
" hit actually being added to TM vector";
838 if (::isDoubleSided(modIdInner, tTopo)) {
843 missingHitAdded =
true;
847 if (innerModInactive) {
848 if (modIdOuter != 0) {
849 LogDebug(
"SiStripHitEfficiency:HitEff") <<
" hit actually being added to TM vector";
852 if (::isDoubleSided(modIdOuter, tTopo)) {
857 missingHitAdded =
true;
864 prev_TKlayers = TKlayers;
869 bool hitsWithBias =
false;
870 for (
auto ilayer : missedLayers) {
871 if (ilayer < TKlayers)
884 const auto nextId = (itm + 1 != TMeas.end()) ? (itm + 1)->recHit()->geographicalId() :
DetId{};
886 if (TKlayers == 9 && theInHit->isValid() && !((!nextId.null()) && (::checkLayer(nextId.rawId(), tTopo) == 9))) {
889 const DetLayer* tob6 = measTracker.geometricSearchTracker()->tobLayers().back();
892 const auto tmp = theLayerMeasurements.measurements(*tob6, tsosTOB5, prop, chi2Estimator);
895 LogDebug(
"SiStripHitEfficiencyWorker") <<
"size of TM from propagation = " <<
tmp.size();
900 const auto& tob6TM =
tmp.back();
901 const auto& tob6Hit = tob6TM.recHit();
902 if (tob6Hit->geographicalId().rawId() != 0) {
903 LogDebug(
"SiStripHitEfficiencyWorker") <<
"tob6 hit actually being added to TM vector";
904 TMs.emplace_back(tob6TM, tTopo, tkgeom,
propagator);
910 if (TKlayers == 21 && theInHit->isValid() &&
911 !((!nextId.null()) && (::checkLayer(nextId.rawId(), tTopo) == 21))) {
912 const DetLayer* tec9pos = measTracker.geometricSearchTracker()->posTecLayers().back();
913 const DetLayer* tec9neg = measTracker.geometricSearchTracker()->negTecLayers().back();
920 LogDebug(
"SiStripHitEfficiencyWorker") <<
"there is a problem with TEC 9 extrapolation";
923 std::vector<TrajectoryMeasurement>
tmp;
924 if (tTopo->tecSide(iidd) == 1) {
925 tmp = theLayerMeasurements.measurements(*tec9neg, tsosTEC9, prop, chi2Estimator);
928 if (tTopo->tecSide(iidd) == 2) {
929 tmp = theLayerMeasurements.measurements(*tec9pos, tsosTEC9, prop, chi2Estimator);
937 const auto& tec9TM =
tmp.back();
938 const auto& tec9Hit = tec9TM.recHit();
940 const unsigned int tec9id = tec9Hit->geographicalId().rawId();
941 LogDebug(
"SiStripHitEfficiencyWorker")
942 <<
"tec9id = " << tec9id <<
" is Double sided = " << ::isDoubleSided(tec9id, tTopo)
943 <<
" and 0x3 = " << (tec9id & 0x3);
945 if (tec9Hit->geographicalId().rawId() != 0) {
946 LogDebug(
"SiStripHitEfficiencyWorker") <<
"tec9 hit actually being added to TM vector";
949 if (::isDoubleSided(tec9id, tTopo)) {
950 TMs.emplace_back(tec9TM, tTopo, tkgeom,
propagator, 1);
951 TMs.emplace_back(tec9TM, tTopo, tkgeom,
propagator, 2);
953 TMs.emplace_back(tec9TM, tTopo, tkgeom,
propagator);
958 for (
const auto& tm : TMs) {
972 LogDebug(
"SiStripHitEfficiencyWorker") <<
"After looping over TrajAtValidHit list";
974 LogDebug(
"SiStripHitEfficiencyWorker") <<
"end TMeasurement loop";
976 LogDebug(
"SiStripHitEfficiencyWorker") <<
"end of trajectories loop";
const edm::ESGetToken< MeasurementTracker, CkfComponentsRecord > measTrackerToken_
const edm::EDGetTokenT< DetIdVector > digisVec_token_
MonitorElement * h_nTracksVsPU
virtual std::vector< DetWithState > compatibleDets(const TrajectoryStateOnSurface &startingState, const Propagator &prop, const MeasurementEstimator &est) const
MonitorElement * h_instLumi
void fillForTraj(const TrajectoryAtInvalidHit &tm, const TrackerTopology *tTopo, const TrackerGeometry *tkgeom, const StripClusterParameterEstimator &stripCPE, const SiStripQuality &stripQuality, const DetIdVector &fedErrorIds, const edm::Handle< edm::DetSetVector< SiStripRawDigi >> &commonModeDigis, const edmNew::DetSetVector< SiStripCluster > &theClusters, int bunchCrossing, float instLumi, float PU, bool highPurity)
const edm::EDGetTokenT< std::vector< Trajectory > > trajectories_token_
std::vector< int > missHitPerLayer
float instLumi() const
Return the luminosity for the current nibble.
std::unordered_map< uint32_t, int > fedErrorCounts
std::unique_ptr< TkHistoMap > FEDErrorOccupancy
const edm::EDGetTokenT< MeasurementTrackerEvent > trackerEvent_token_
const edm::ESGetToken< StripClusterParameterEstimator, TkStripCPERecord > stripCPEToken_
const edm::ESGetToken< Propagator, TrackingComponentsRecord > propagatorToken_
const edm::EDGetTokenT< OnlineLuminosityRecord > metaDataToken_
const edm::ESGetToken< Chi2MeasurementEstimatorBase, TrackingComponentsRecord > chi2EstimatorToken_
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
unsigned int trackMultiplicityCut_
bool useAllHitsFromTracksWithMissingHits_
dqm::reco::MonitorElement * EventStats
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
const edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > magFieldToken_
SiStripHitEffData calibData_
const edm::EDGetTokenT< TrajTrackAssociationCollection > trajTrackAsso_token_
const edm::EDGetTokenT< DetIdCollection > digisCol_token_
MonitorElement * h_nTracks
const edm::EDGetTokenT< edmNew::DetSetVector< SiStripCluster > > clusters_token_
const edm::EDGetTokenT< reco::TrackCollection > combinatorialTracks_token_
const edm::ESGetToken< SiStripQuality, SiStripQualityRcd > stripQualityToken_
const edm::EDGetTokenT< LumiScalersCollection > scalerToken_
std::vector< unsigned int > hitRecoveryCounters
const edm::EDGetTokenT< edm::DetSetVector< SiStripRawDigi > > commonModeToken_
Log< level::Warning, false > LogWarning
float avgPileUp() const
Return the average pileup for th current nibble.
ConstRecHitPointer const & recHit() const
bool doMissingHitsRecovery_