Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:59:55

0001 ////////////////////////////////////////////////////////////////////////////////
0002 // Package:          CalibTracker/SiStripHitEfficiency
0003 // Class:            SiStripHitEfficiencyWorker
0004 // Original Author:  Pieter David
0005 //
0006 // Adapted from      HitEff (Keith Ulmer -- University of Colorado, keith.ulmer@colorado.edu
0007 //                   SiStripHitEffFromCalibTree (Christopher Edelmaier)
0008 ///////////////////////////////////////////////////////////////////////////////
0009 
0010 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0011 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0012 #include "CalibTracker/SiStripHitEfficiency/interface/SiStripHitEffData.h"
0013 #include "CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h"
0014 #include "CalibTracker/SiStripHitEfficiency/interface/TrajectoryAtInvalidHit.h"
0015 #include "DQM/SiStripCommon/interface/TkHistoMap.h"
0016 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0017 #include "DataFormats/Common/interface/DetSetVector.h"
0018 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0019 #include "DataFormats/Common/interface/Handle.h"
0020 #include "DataFormats/DetId/interface/DetIdCollection.h"
0021 #include "DataFormats/DetId/interface/DetIdVector.h"
0022 #include "DataFormats/GeometryCommonDetAlgo/interface/MeasurementError.h"
0023 #include "DataFormats/GeometryCommonDetAlgo/interface/MeasurementVector.h"
0024 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
0025 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0026 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0027 #include "DataFormats/GeometryVector/interface/LocalVector.h"
0028 #include "DataFormats/OnlineMetaData/interface/OnlineLuminosityRecord.h"
0029 #include "DataFormats/Scalers/interface/LumiScalers.h"
0030 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0031 #include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" /* for STRIPS_PER_APV*/
0032 #include "DataFormats/SiStripDigi/interface/SiStripRawDigi.h"
0033 #include "DataFormats/TrackReco/interface/Track.h"
0034 #include "DataFormats/TrackReco/interface/TrackBase.h"
0035 #include "DataFormats/TrackReco/interface/TrackExtra.h"
0036 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0037 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0038 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0039 #include "FWCore/Framework/interface/Event.h"
0040 #include "FWCore/Framework/interface/EventSetup.h"
0041 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0042 #include "FWCore/ParameterSet/interface/ParameterDescription.h"
0043 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0044 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0045 #include "FWCore/Utilities/interface/Exception.h"
0046 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0047 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0048 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0049 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0050 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0051 #include "MagneticField/Engine/interface/MagneticField.h"
0052 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0053 #include "RecoTracker/MeasurementDet/interface/MeasurementTracker.h"
0054 #include "RecoTracker/MeasurementDet/interface/MeasurementTrackerEvent.h"
0055 #include "RecoTracker/Record/interface/CkfComponentsRecord.h"
0056 #include "TrackingTools/DetLayers/interface/DetLayer.h"
0057 #include "TrackingTools/GeomPropagators/interface/AnalyticalPropagator.h"
0058 #include "TrackingTools/KalmanUpdators/interface/Chi2MeasurementEstimator.h"
0059 #include "TrackingTools/MaterialEffects/interface/PropagatorWithMaterial.h"
0060 #include "TrackingTools/MeasurementDet/interface/LayerMeasurements.h"
0061 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0062 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0063 
0064 class SiStripHitEfficiencyWorker : public DQMEDAnalyzer {
0065 public:
0066   explicit SiStripHitEfficiencyWorker(const edm::ParameterSet& conf);
0067   ~SiStripHitEfficiencyWorker() override = default;
0068   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0069 
0070 private:
0071   void bookHistograms(DQMStore::IBooker& booker, const edm::Run& run, const edm::EventSetup& setup) override;
0072   void analyze(const edm::Event& e, const edm::EventSetup& c) override;
0073   void fillForTraj(const TrajectoryAtInvalidHit& tm,
0074                    const TrackerTopology* tTopo,
0075                    const TrackerGeometry* tkgeom,
0076                    const StripClusterParameterEstimator& stripCPE,
0077                    const SiStripQuality& stripQuality,
0078                    const DetIdVector& fedErrorIds,
0079                    const edm::Handle<edm::DetSetVector<SiStripRawDigi>>& commonModeDigis,
0080                    const edmNew::DetSetVector<SiStripCluster>& theClusters,
0081                    int bunchCrossing,
0082                    float instLumi,
0083                    float PU,
0084                    bool highPurity);
0085 
0086   // ----------member data ---------------------------
0087   SiStripHitEffData calibData_;
0088 
0089   // event data tokens
0090   const edm::EDGetTokenT<LumiScalersCollection> scalerToken_;
0091   const edm::EDGetTokenT<OnlineLuminosityRecord> metaDataToken_;
0092   const edm::EDGetTokenT<edm::DetSetVector<SiStripRawDigi>> commonModeToken_;
0093   const edm::EDGetTokenT<reco::TrackCollection> combinatorialTracks_token_;
0094   const edm::EDGetTokenT<std::vector<Trajectory>> trajectories_token_;
0095   const edm::EDGetTokenT<TrajTrackAssociationCollection> trajTrackAsso_token_;
0096   const edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster>> clusters_token_;
0097   const edm::EDGetTokenT<DetIdCollection> digisCol_token_;
0098   const edm::EDGetTokenT<DetIdVector> digisVec_token_;
0099   const edm::EDGetTokenT<MeasurementTrackerEvent> trackerEvent_token_;
0100 
0101   // event setup tokens
0102   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0103   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_;
0104   const edm::ESGetToken<StripClusterParameterEstimator, TkStripCPERecord> stripCPEToken_;
0105   const edm::ESGetToken<SiStripQuality, SiStripQualityRcd> stripQualityToken_;
0106   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magFieldToken_;
0107   const edm::ESGetToken<MeasurementTracker, CkfComponentsRecord> measTrackerToken_;
0108   const edm::ESGetToken<Chi2MeasurementEstimatorBase, TrackingComponentsRecord> chi2EstimatorToken_;
0109   const edm::ESGetToken<Propagator, TrackingComponentsRecord> propagatorToken_;
0110   const edm::ESGetToken<TkDetMap, TrackerTopologyRcd> tkDetMapToken_;
0111 
0112   // configurable parameters
0113   std::string dqmDir_;
0114   unsigned int layers_;
0115   bool DEBUG_;
0116   bool addLumi_;
0117   bool addCommonMode_;
0118   bool cutOnTracks_;
0119   unsigned int trackMultiplicityCut_;
0120   bool useFirstMeas_;
0121   bool useLastMeas_;
0122   bool useAllHitsFromTracksWithMissingHits_;
0123   bool doMissingHitsRecovery_;
0124   unsigned int clusterMatchingMethod_;
0125   float resXSig_;
0126   float clusterTracjDist_;
0127   float stripsApvEdge_;
0128   bool useOnlyHighPurityTracks_;
0129   int bunchX_;
0130   bool showRings_;
0131   bool showTOB6TEC9_;
0132   unsigned int nTEClayers_;
0133 
0134   // output file
0135   std::set<uint32_t> badModules_;
0136 
0137   // for the missing hits recovery
0138   std::vector<unsigned int> hitRecoveryCounters;
0139   std::vector<unsigned int> hitTotalCounters;
0140   int totalNbHits;
0141   std::vector<int> missHitPerLayer;
0142 
0143   struct EffME1 {
0144     EffME1() : hTotal(nullptr), hFound(nullptr) {}
0145     EffME1(MonitorElement* total, MonitorElement* found) : hTotal(total), hFound(found) {}
0146 
0147     void fill(double x, bool found, float weight = 1.) {
0148       hTotal->Fill(x, weight);
0149       if (found) {
0150         hFound->Fill(x, weight);
0151       }
0152     }
0153 
0154     MonitorElement *hTotal, *hFound;
0155   };
0156   struct EffTkMap {
0157     EffTkMap() : hTotal(nullptr), hFound(nullptr) {}
0158     EffTkMap(std::unique_ptr<TkHistoMap>&& total, std::unique_ptr<TkHistoMap>&& found)
0159         : hTotal(std::move(total)), hFound(std::move(found)) {}
0160 
0161     void fill(uint32_t id, bool found, float weight = 1.) {
0162       hTotal->fill(id, weight);
0163       if (found) {
0164         hFound->fill(id, weight);
0165       }
0166     }
0167 
0168     bool check(uint32_t id) {
0169       if (hTotal->getValue(id) < hFound->getValue(id)) {
0170         return false;
0171       } else {
0172         return true;
0173       }
0174     }
0175 
0176     std::unique_ptr<TkHistoMap> hTotal, hFound;
0177   };
0178 
0179   MonitorElement *h_bx, *h_instLumi, *h_PU;
0180   MonitorElement *h_nTracks, *h_nTracksVsPU;
0181   EffME1 h_goodLayer;
0182   EffME1 h_allLayer;
0183   EffME1 h_layer;
0184   std::vector<MonitorElement*> h_resolution;
0185   std::vector<EffME1> h_layer_vsLumi;
0186   std::vector<EffME1> h_layer_vsBx;
0187   std::vector<EffME1> h_layer_vsPU;
0188   std::vector<EffME1> h_layer_vsCM;
0189   std::vector<MonitorElement*> h_hotcold;
0190 
0191   EffTkMap h_module;
0192 };
0193 
0194 //
0195 // constructors and destructor
0196 //
0197 
0198 SiStripHitEfficiencyWorker::SiStripHitEfficiencyWorker(const edm::ParameterSet& conf)
0199     : scalerToken_(consumes<LumiScalersCollection>(conf.getParameter<edm::InputTag>("lumiScalers"))),
0200       metaDataToken_(consumes<OnlineLuminosityRecord>(conf.getParameter<edm::InputTag>("metadata"))),
0201       commonModeToken_(mayConsume<edm::DetSetVector<SiStripRawDigi>>(conf.getParameter<edm::InputTag>("commonMode"))),
0202       combinatorialTracks_token_(
0203           consumes<reco::TrackCollection>(conf.getParameter<edm::InputTag>("combinatorialTracks"))),
0204       trajectories_token_(consumes<std::vector<Trajectory>>(conf.getParameter<edm::InputTag>("trajectories"))),
0205       trajTrackAsso_token_(consumes<TrajTrackAssociationCollection>(conf.getParameter<edm::InputTag>("trajectories"))),
0206       clusters_token_(
0207           consumes<edmNew::DetSetVector<SiStripCluster>>(conf.getParameter<edm::InputTag>("siStripClusters"))),
0208       digisCol_token_(consumes(conf.getParameter<edm::InputTag>("siStripDigis"))),
0209       digisVec_token_(consumes(conf.getParameter<edm::InputTag>("siStripDigis"))),
0210       trackerEvent_token_(consumes<MeasurementTrackerEvent>(conf.getParameter<edm::InputTag>("trackerEvent"))),
0211       tTopoToken_(esConsumes()),
0212       tkGeomToken_(esConsumes()),
0213       stripCPEToken_(esConsumes(edm::ESInputTag{"", "StripCPEfromTrackAngle"})),
0214       stripQualityToken_(esConsumes()),
0215       magFieldToken_(esConsumes()),
0216       measTrackerToken_(esConsumes()),
0217       chi2EstimatorToken_(esConsumes(edm::ESInputTag{"", "Chi2"})),
0218       propagatorToken_(esConsumes(edm::ESInputTag{"", "PropagatorWithMaterial"})),
0219       tkDetMapToken_(esConsumes<edm::Transition::BeginRun>()),
0220       dqmDir_(conf.getParameter<std::string>("dqmDir")),
0221       layers_(conf.getParameter<int>("Layer")),
0222       DEBUG_(conf.getUntrackedParameter<bool>("Debug", false)),
0223       addLumi_(conf.getUntrackedParameter<bool>("addLumi", false)),
0224       addCommonMode_(conf.getUntrackedParameter<bool>("addCommonMode", false)),
0225       cutOnTracks_(conf.getParameter<bool>("cutOnTracks")),
0226       trackMultiplicityCut_(conf.getParameter<unsigned int>("trackMultiplicity")),
0227       useFirstMeas_(conf.getParameter<bool>("useFirstMeas")),
0228       useLastMeas_(conf.getParameter<bool>("useLastMeas")),
0229       useAllHitsFromTracksWithMissingHits_(conf.getParameter<bool>("useAllHitsFromTracksWithMissingHits")),
0230       doMissingHitsRecovery_(conf.getParameter<bool>("doMissingHitsRecovery")),
0231       clusterMatchingMethod_(conf.getParameter<int>("ClusterMatchingMethod")),
0232       resXSig_(conf.getParameter<double>("ResXSig")),
0233       clusterTracjDist_(conf.getParameter<double>("ClusterTrajDist")),
0234       stripsApvEdge_(conf.getParameter<double>("StripsApvEdge")),
0235       useOnlyHighPurityTracks_(conf.getParameter<bool>("UseOnlyHighPurityTracks")),
0236       bunchX_(conf.getUntrackedParameter<int>("BunchCrossing", 0)),
0237       showRings_(conf.getUntrackedParameter<bool>("ShowRings", false)),
0238       showTOB6TEC9_(conf.getUntrackedParameter<bool>("ShowTOB6TEC9", false)) {
0239   hitRecoveryCounters.resize(k_END_OF_LAYERS, 0);
0240   hitTotalCounters.resize(k_END_OF_LAYERS, 0);
0241   missHitPerLayer.resize(k_END_OF_LAYERS, 0);
0242   totalNbHits = 0;
0243 
0244   nTEClayers_ = (showRings_ ? 7 : 9);  // number of rings or wheels
0245 
0246   const std::string badModulesFile = conf.getUntrackedParameter<std::string>("BadModulesFile", "");
0247   if (!badModulesFile.empty()) {
0248     std::ifstream badModules_file(badModulesFile);
0249     uint32_t badmodule_detid;
0250     int mods, fiber1, fiber2, fiber3;
0251     if (badModules_file.is_open()) {
0252       std::string line;
0253       while (getline(badModules_file, line)) {
0254         if (badModules_file.eof())
0255           continue;
0256         std::stringstream ss(line);
0257         ss >> badmodule_detid >> mods >> fiber1 >> fiber2 >> fiber3;
0258         if (badmodule_detid != 0 && mods == 1 && (fiber1 == 1 || fiber2 == 1 || fiber3 == 1))
0259           badModules_.insert(badmodule_detid);
0260       }
0261       badModules_file.close();
0262     }
0263   }
0264   if (!badModules_.empty())
0265     LogDebug("SiStripHitEfficiencyWorker") << "Remove additionnal bad modules from the analysis: ";
0266   for (const auto badMod : badModules_) {
0267     LogDebug("SiStripHitEfficiencyWorker") << " " << badMod;
0268   }
0269 }
0270 
0271 void SiStripHitEfficiencyWorker::bookHistograms(DQMStore::IBooker& booker,
0272                                                 const edm::Run& run,
0273                                                 const edm::EventSetup& setup) {
0274   booker.setCurrentFolder(fmt::format("{}/EventInfo", dqmDir_));
0275   h_bx = booker.book1D("bx", "bx", 3600, 0, 3600);
0276   h_instLumi = booker.book1D("instLumi", "inst. lumi.", 250, 0, 25000);
0277   h_PU = booker.book1D("PU", "PU", 200, 0, 200);
0278   h_nTracks = booker.book1D("ntracks", "n.tracks;n. tracks;n.events", 500, -0.5, 499.5);
0279   h_nTracksVsPU = booker.bookProfile("nTracksVsPU", "n. tracks vs PU; PU; n.tracks ", 200, 0, 200, 500, -0.5, 499.5);
0280 
0281   calibData_.EventStats = booker.book2I("EventStats", "Statistics", 3, -0.5, 2.5, 1, 0, 1);
0282   calibData_.EventStats->setBinLabel(1, "events count", 1);
0283   calibData_.EventStats->setBinLabel(2, "tracks count", 1);
0284   calibData_.EventStats->setBinLabel(3, "measurements count", 1);
0285 
0286   booker.setCurrentFolder(dqmDir_);
0287   h_goodLayer = EffME1(booker.book1D("goodlayer_total", "goodlayer_total", 35, 0., 35.),
0288                        booker.book1D("goodlayer_found", "goodlayer_found", 35, 0., 35.));
0289   h_allLayer = EffME1(booker.book1D("alllayer_total", "alllayer_total", 35, 0., 35.),
0290                       booker.book1D("alllayer_found", "alllayer_found", 35, 0., 35.));
0291 
0292   h_layer = EffME1(
0293       booker.book1D(
0294           "layer_found", "layer_found", bounds::k_END_OF_LAYERS, 0., static_cast<float>(bounds::k_END_OF_LAYERS)),
0295       booker.book1D(
0296           "layer_total", "layer_total", bounds::k_END_OF_LAYERS, 0., static_cast<float>(bounds::k_END_OF_LAYERS)));
0297 
0298   for (int layer = 1; layer != bounds::k_END_OF_LAYERS; ++layer) {
0299     const auto lyrName = ::layerName(layer, showRings_, nTEClayers_);
0300 
0301     // book resolutions
0302     booker.setCurrentFolder(fmt::format("{}/Resolutions", dqmDir_));
0303     auto ihres = booker.book1D(Form("resol_layer_%i", layer), lyrName, 125, -125., 125.);
0304     ihres->setAxisTitle("trajX-clusX [strip unit]");
0305     h_resolution.push_back(ihres);
0306 
0307     // book plots vs Lumi
0308     booker.setCurrentFolder(fmt::format("{}/VsLumi", dqmDir_));
0309     h_layer_vsLumi.push_back(EffME1(booker.book1D(Form("layertotal_vsLumi_layer_%i", layer), lyrName, 100, 0, 25000),
0310                                     booker.book1D(Form("layerfound_vsLumi_layer_%i", layer), lyrName, 100, 0, 25000)));
0311 
0312     // book plots vs Lumi
0313     booker.setCurrentFolder(fmt::format("{}/VsPu", dqmDir_));
0314     h_layer_vsPU.push_back(EffME1(booker.book1D(Form("layertotal_vsPU_layer_%i", layer), lyrName, 45, 0, 90),
0315                                   booker.book1D(Form("layerfound_vsPU_layer_%i", layer), lyrName, 45, 0, 90)));
0316     if (addCommonMode_) {
0317       // book plots for common mode
0318       booker.setCurrentFolder(fmt::format("{}/CommonMode", dqmDir_));
0319       h_layer_vsCM.push_back(EffME1(booker.book1D(Form("layertotal_vsCM_layer_%i", layer), lyrName, 20, 0, 400),
0320                                     booker.book1D(Form("layerfound_vsCM_layer_%i", layer), lyrName, 20, 0, 400)));
0321     }
0322 
0323     // book plots vs Lumi
0324     booker.setCurrentFolder(fmt::format("{}/VsBx", dqmDir_));
0325     h_layer_vsBx.push_back(EffME1(
0326         booker.book1D(Form("totalVsBx_layer%i", layer), Form("layer %i (%s)", layer, lyrName.c_str()), 3565, 0, 3565),
0327         booker.book1D(Form("foundVsBx_layer%i", layer), Form("layer %i (%s)", layer, lyrName.c_str()), 3565, 0, 3565)));
0328 
0329     // book hot and cold
0330     booker.setCurrentFolder(fmt::format("{}/MissingHits", dqmDir_));
0331     if (layer <= bounds::k_LayersAtTOBEnd) {
0332       const bool isTIB = layer <= bounds::k_LayersAtTIBEnd;
0333       const auto partition = (isTIB ? "TIB" : "TOB");
0334       const auto yMax = (isTIB ? 100 : 120);
0335 
0336       const auto& tit =
0337           Form("%s%i: Map of missing hits", partition, (isTIB ? layer : layer - bounds::k_LayersAtTIBEnd));
0338 
0339       // histogram name must not contain ":" otherwise it fails upload to the GUI
0340       // see https://github.com/cms-DQM/dqmgui_prod/blob/af0a388e8f57c60e51111585d298aeeea943367f/src/cpp/DQM/DQMStore.cc#L56
0341       std::string name{tit};
0342       ::replaceInString(name, ":", "");
0343 
0344       auto ihhotcold = booker.book2D(name, tit, 100, -1, 361, 100, -yMax, yMax);
0345       ihhotcold->setAxisTitle("#phi [deg]", 1);
0346       ihhotcold->setBinLabel(1, "360", 1);
0347       ihhotcold->setBinLabel(50, "180", 1);
0348       ihhotcold->setBinLabel(100, "0", 1);
0349       ihhotcold->setAxisTitle("Global Z [cm]", 2);
0350       ihhotcold->setOption("colz");
0351       h_hotcold.push_back(ihhotcold);
0352     } else {
0353       const bool isTID = layer <= bounds::k_LayersAtTIDEnd;
0354       const auto partitions =
0355           (isTID ? std::vector<std::string>{"TIDplus", "TIDminus"} : std::vector<std::string>{"TECplus", "TECminus"});
0356       const auto axMax = (isTID ? 100 : 120);
0357       for (const auto& part : partitions) {
0358         // create the title by replacing the minus/plus symbols
0359         std::string forTitle{part};
0360         ::replaceInString(forTitle, "minus", "-");
0361         ::replaceInString(forTitle, "plus", "+");
0362 
0363         // histogram name must not contain ":" otherwise it fails upload to the GUI
0364         // see https://github.com/cms-DQM/dqmgui_prod/blob/af0a388e8f57c60e51111585d298aeeea943367f/src/cpp/DQM/DQMStore.cc#L56
0365         const auto& name = Form("%s %i Map of Missing Hits",
0366                                 part.c_str(),
0367                                 (isTID ? layer - bounds::k_LayersAtTOBEnd : layer - bounds::k_LayersAtTIDEnd));
0368         const auto& tit = Form("%s%i: Map of Missing Hits",
0369                                forTitle.c_str(),
0370                                (isTID ? layer - bounds::k_LayersAtTOBEnd : layer - bounds::k_LayersAtTIDEnd));
0371 
0372         auto ihhotcold = booker.book2D(name, tit, 100, -axMax, axMax, 100, -axMax, axMax);
0373         ihhotcold->setAxisTitle("Global Y", 1);
0374         ihhotcold->setBinLabel(1, "+Y", 1);
0375         ihhotcold->setBinLabel(50, "0", 1);
0376         ihhotcold->setBinLabel(100, "-Y", 1);
0377         ihhotcold->setAxisTitle("Global X", 2);
0378         ihhotcold->setBinLabel(1, "-X", 2);
0379         ihhotcold->setBinLabel(50, "0", 2);
0380         ihhotcold->setBinLabel(100, "+X", 2);
0381         ihhotcold->setOption("colz");
0382         h_hotcold.push_back(ihhotcold);
0383       }
0384     }
0385   }
0386 
0387   // come back to the main folder
0388   booker.setCurrentFolder(dqmDir_);
0389   const auto tkDetMapFolder = fmt::format("{}/TkDetMaps", dqmDir_);
0390 
0391   const TkDetMap* tkDetMap = &setup.getData(tkDetMapToken_);
0392   h_module =
0393       EffTkMap(std::make_unique<TkHistoMap>(tkDetMap, booker, tkDetMapFolder, "perModule_total", 0, false, true),
0394                std::make_unique<TkHistoMap>(tkDetMap, booker, tkDetMapFolder, "perModule_found", 0, false, true));
0395 
0396   // fill the FED Errors
0397   booker.setCurrentFolder(dqmDir_);
0398   const auto FEDErrorMapFolder = fmt::format("{}/FEDErrorTkDetMaps", dqmDir_);
0399   calibData_.FEDErrorOccupancy =
0400       std::make_unique<TkHistoMap>(tkDetMap, booker, FEDErrorMapFolder, "perModule_FEDErrors", 0, false, true);
0401 }
0402 
0403 void SiStripHitEfficiencyWorker::analyze(const edm::Event& e, const edm::EventSetup& es) {
0404   const auto tTopo = &es.getData(tTopoToken_);
0405 
0406   //  bool DEBUG_ = false;
0407 
0408   LogDebug("SiStripHitEfficiencyWorker") << "beginning analyze from HitEff";
0409 
0410   // Step A: Get Inputs
0411 
0412   // Luminosity informations
0413   edm::Handle<LumiScalersCollection> lumiScalers = e.getHandle(scalerToken_);
0414   edm::Handle<OnlineLuminosityRecord> metaData = e.getHandle(metaDataToken_);
0415 
0416   float instLumi = 0;
0417   float PU = 0;
0418   if (addLumi_) {
0419     if (lumiScalers.isValid() && !lumiScalers->empty()) {
0420       if (lumiScalers->begin() != lumiScalers->end()) {
0421         instLumi = lumiScalers->begin()->instantLumi();
0422         PU = lumiScalers->begin()->pileup();
0423       }
0424     } else if (metaData.isValid()) {
0425       instLumi = metaData->instLumi();
0426       PU = metaData->avgPileUp();
0427     } else {
0428       edm::LogWarning("SiStripHitEfficiencyWorker") << "could not find a source for the Luminosity and PU";
0429     }
0430   }
0431 
0432   h_bx->Fill(e.bunchCrossing());
0433   h_instLumi->Fill(instLumi);
0434   h_PU->Fill(PU);
0435 
0436   edm::Handle<edm::DetSetVector<SiStripRawDigi>> commonModeDigis;
0437   if (addCommonMode_)
0438     e.getByToken(commonModeToken_, commonModeDigis);
0439 
0440   edm::Handle<reco::TrackCollection> tracksCKF;
0441   e.getByToken(combinatorialTracks_token_, tracksCKF);
0442 
0443   edm::Handle<std::vector<Trajectory>> TrajectoryCollectionCKF;
0444   e.getByToken(trajectories_token_, TrajectoryCollectionCKF);
0445 
0446   edm::Handle<TrajTrackAssociationCollection> trajTrackAssociationHandle;
0447   e.getByToken(trajTrackAsso_token_, trajTrackAssociationHandle);
0448 
0449   edm::Handle<edmNew::DetSetVector<SiStripCluster>> theClusters;
0450   e.getByToken(clusters_token_, theClusters);
0451 
0452   // get the list of module IDs with FED-detected errors
0453   //  - In Aug-2023, the data format was changed from DetIdCollection to DetIdVector.
0454   //  - To provide some level of backward-compatibility,
0455   //    the plugin checks for both types giving preference to the new format.
0456   //  - If only the old format is available, the collection is
0457   //    converted to the new format, then used downstream.
0458   auto const& fedErrorIdsCol_h = e.getHandle(digisCol_token_);
0459   auto const& fedErrorIdsVec_h = e.getHandle(digisVec_token_);
0460   if (not fedErrorIdsCol_h.isValid() and not fedErrorIdsVec_h.isValid()) {
0461     throw cms::Exception("InvalidProductSiStripDetIdsWithFEDErrors")
0462         << "no valid product for SiStrip DetIds with FED errors (see parameter \"siStripDigis\"), "
0463            "neither for new format (DetIdVector) nor old format (DetIdCollection)";
0464   }
0465   auto const& fedErrorIds = fedErrorIdsVec_h.isValid() ? *fedErrorIdsVec_h : fedErrorIdsCol_h->as_vector();
0466 
0467   // fill the calibData with the FEDErrors
0468   for (const auto& fedErr : fedErrorIds) {
0469     // fill the TkHistoMap occupancy map
0470     calibData_.FEDErrorOccupancy->fill(fedErr.rawId(), 1.);
0471 
0472     // fill the unordered map
0473     if (calibData_.fedErrorCounts.find(fedErr.rawId()) != calibData_.fedErrorCounts.end()) {
0474       calibData_.fedErrorCounts[fedErr.rawId()] += 1;
0475     } else {
0476       calibData_.fedErrorCounts.insert(std::make_pair(fedErr.rawId(), 1));
0477     }
0478   }
0479 
0480   edm::Handle<MeasurementTrackerEvent> measurementTrackerEvent;
0481   e.getByToken(trackerEvent_token_, measurementTrackerEvent);
0482 
0483   const auto tkgeom = &es.getData(tkGeomToken_);
0484   const auto& stripcpe = es.getData(stripCPEToken_);
0485   const auto& stripQuality = es.getData(stripQualityToken_);
0486   const auto& magField = es.getData(magFieldToken_);
0487   const auto& measTracker = es.getData(measTrackerToken_);
0488   const auto& chi2Estimator = es.getData(chi2EstimatorToken_);
0489   const auto& prop = es.getData(propagatorToken_);
0490 
0491   // Tracking
0492   LogDebug("SiStripHitEfficiencyWorker") << "number ckf tracks found = " << tracksCKF->size();
0493 
0494   h_nTracks->Fill(tracksCKF->size());
0495   h_nTracksVsPU->Fill(PU, tracksCKF->size());
0496 
0497   // bin 0: one entry for each event
0498   calibData_.EventStats->Fill(0., 0., 1);
0499   // bin 1: one entry for each track
0500   calibData_.EventStats->Fill(1., 0., tracksCKF->size());
0501 
0502   if (!tracksCKF->empty()) {
0503     if (cutOnTracks_ && (tracksCKF->size() >= trackMultiplicityCut_))
0504       return;
0505     if (cutOnTracks_)
0506       LogDebug("SiStripHitEfficiencyWorker")
0507           << "starting checking good event with < " << trackMultiplicityCut_ << " tracks";
0508 
0509     // actually should do a loop over all the tracks in the event here
0510 
0511     // Looping over traj-track associations to be able to get traj & track informations
0512     for (const auto& trajTrack : *trajTrackAssociationHandle) {
0513       // for each track, fill some variables such as number of hits and momentum
0514 
0515       const bool highPurity = trajTrack.val->quality(reco::TrackBase::TrackQuality::highPurity);
0516       auto TMeas = trajTrack.key->measurements();
0517       totalNbHits += int(TMeas.size());
0518 
0519       /*
0520       const bool hasMissingHits = std::any_of(std::begin(TMeas), std::end(TMeas), [](const auto& tm) {
0521         return tm.recHit()->getType() == TrackingRecHit::Type::missing;
0522       });
0523       */
0524 
0525       // Check whether the trajectory has some missing hits
0526       bool hasMissingHits{false};
0527       int previous_layer{999};
0528       std::vector<unsigned int> missedLayers;
0529 
0530       for (const auto& itm : TMeas) {
0531         auto theHit = itm.recHit();
0532         unsigned int iidd = theHit->geographicalId().rawId();
0533         int layer = ::checkLayer(iidd, tTopo);
0534         int missedLayer = layer + 1;
0535         int diffPreviousLayer = (layer - previous_layer);
0536         if (doMissingHitsRecovery_) {
0537           //Layers from TIB + TOB
0538           if (diffPreviousLayer == -2 && missedLayer > k_LayersStart && missedLayer < k_LayersAtTOBEnd) {
0539             missHitPerLayer[missedLayer] += 1;
0540             hasMissingHits = true;
0541           }
0542           //Layers from TID
0543           else if (diffPreviousLayer == -2 && (missedLayer > k_LayersAtTOBEnd + 1 && missedLayer <= k_LayersAtTIDEnd)) {
0544             missHitPerLayer[missedLayer] += 1;
0545             hasMissingHits = true;
0546           }
0547           //Layers from TEC
0548           else if (diffPreviousLayer == -2 && missedLayer > k_LayersAtTIDEnd && missedLayer <= k_LayersAtTECEnd) {
0549             missHitPerLayer[missedLayer] += 1;
0550             hasMissingHits = true;
0551           }
0552 
0553           //##### TID Layer 11 in transition TID -> TIB : layer is in TIB, previous layer  = 12
0554           if ((layer > k_LayersStart && layer <= k_LayersAtTIBEnd) && (previous_layer == 12)) {
0555             missHitPerLayer[11] += 1;
0556             hasMissingHits = true;
0557           }
0558 
0559           //##### TEC Layer 14 in transition TEC -> TOB : layer is in TOB, previous layer =  15
0560           if ((layer > k_LayersAtTIBEnd && layer <= k_LayersAtTOBEnd) && (previous_layer == 15)) {
0561             missHitPerLayer[14] += 1;
0562             hasMissingHits = true;
0563           }
0564         }
0565         if (theHit->getType() == TrackingRecHit::Type::missing)
0566           hasMissingHits = true;
0567 
0568         if (hasMissingHits)
0569           missedLayers.push_back(layer);
0570         previous_layer = layer;
0571       }
0572 
0573       // Loop on each measurement and take into consideration
0574       //--------------------------------------------------------
0575       unsigned int prev_TKlayers = 0;
0576       for (auto itm = TMeas.cbegin(); itm != TMeas.cend(); ++itm) {
0577         const auto theInHit = (*itm).recHit();
0578 
0579         //bin 2: one entry for each measurement
0580         calibData_.EventStats->Fill(2., 0., 1.);
0581 
0582         LogDebug("SiStripHitEfficiencyWorker") << "theInHit is valid = " << theInHit->isValid();
0583 
0584         unsigned int iidd = theInHit->geographicalId().rawId();
0585 
0586         unsigned int TKlayers = ::checkLayer(iidd, tTopo);
0587 
0588         // do not bother with pixel hits
0589         if (DetId(iidd).subdetId() < SiStripSubdetector::TIB)
0590           continue;
0591 
0592         LogDebug("SiStripHitEfficiencyWorker") << "TKlayer from trajectory: " << TKlayers << "  from module = " << iidd
0593                                                << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/"
0594                                                << ((iidd & 0x3) == 1) << "/" << ((iidd & 0x3) == 2);
0595 
0596         // Test first and last points of the trajectory
0597         // the list of measurements starts from outer layers  !!! This could change -> should add a check
0598         bool isFirstMeas = (itm == (TMeas.end() - 1));
0599         bool isLastMeas = (itm == (TMeas.begin()));
0600 
0601         if (!useFirstMeas_ && isFirstMeas)
0602           continue;
0603         if (!useLastMeas_ && isLastMeas)
0604           continue;
0605 
0606         // In case of missing hit in the track, check whether to use the other hits or not.
0607         if (hasMissingHits && theInHit->getType() != TrackingRecHit::Type::missing &&
0608             !useAllHitsFromTracksWithMissingHits_)
0609           continue;
0610 
0611         // If Trajectory measurement from TOB 6 or TEC 9, skip it because it's always valid they are filled later
0612         if (TKlayers == bounds::k_LayersAtTOBEnd || TKlayers == bounds::k_LayersAtTECEnd) {
0613           LogDebug("SiStripHitEfficiencyWorker") << "skipping original TM for TOB 6 or TEC 9";
0614           continue;
0615         }
0616 
0617         std::vector<TrajectoryAtInvalidHit> TMs;
0618 
0619         // Make AnalyticalPropagat // TODO where to save these?or to use in TAVH constructor
0620         AnalyticalPropagator propagator(&magField, anyDirection);
0621 
0622         // for double sided layers check both sensors--if no hit was found on either sensor surface,
0623         // the trajectory measurements only have one invalid hit entry on the matched surface
0624         // so get the TrajectoryAtInvalidHit for both surfaces and include them in the study
0625         if (::isDoubleSided(iidd, tTopo) && ((iidd & 0x3) == 0)) {
0626           // do hit eff check twice--once for each sensor
0627           //add a TM for each surface
0628           TMs.emplace_back(*itm, tTopo, tkgeom, propagator, 1);
0629           TMs.emplace_back(*itm, tTopo, tkgeom, propagator, 2);
0630         } else if (::isDoubleSided(iidd, tTopo) && (!::check2DPartner(iidd, TMeas))) {
0631           // if only one hit was found the trajectory measurement is on that sensor surface, and the other surface from
0632           // the matched layer should be added to the study as well
0633           TMs.emplace_back(*itm, tTopo, tkgeom, propagator, 1);
0634           TMs.emplace_back(*itm, tTopo, tkgeom, propagator, 2);
0635           LogDebug("SiStripHitEfficiencyWorker") << " found a hit with a missing partner";
0636         } else {
0637           //only add one TM for the single surface and the other will be added in the next iteration
0638           TMs.emplace_back(*itm, tTopo, tkgeom, propagator);
0639         }
0640 
0641         bool missingHitAdded{false};
0642         std::vector<TrajectoryMeasurement> tmpTmeas;
0643         unsigned int misLayer = TKlayers + 1;
0644         //Use bool doMissingHitsRecovery to add possible missing hits based on actual/previous hit
0645         if (doMissingHitsRecovery_) {
0646           if (int(TKlayers) - int(prev_TKlayers) == -2) {
0647             const DetLayer* detlayer = itm->layer();
0648             const LayerMeasurements layerMeasurements{measTracker, *measurementTrackerEvent};
0649             const TrajectoryStateOnSurface tsos = itm->updatedState();
0650             std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, prop, chi2Estimator);
0651 
0652             if (misLayer > k_LayersAtTIDEnd && misLayer < k_LayersAtTECEnd) {  //TEC
0653               std::vector<ForwardDetLayer const*> negTECLayers = measTracker.geometricSearchTracker()->negTecLayers();
0654               std::vector<ForwardDetLayer const*> posTECLayers = measTracker.geometricSearchTracker()->posTecLayers();
0655               const DetLayer* tecLayerneg = negTECLayers[misLayer - k_LayersAtTIDEnd - 1];
0656               const DetLayer* tecLayerpos = posTECLayers[misLayer - k_LayersAtTIDEnd - 1];
0657               if (tTopo->tecSide(iidd) == 1) {
0658                 tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, prop, chi2Estimator);
0659               } else if (tTopo->tecSide(iidd) == 2) {
0660                 tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, prop, chi2Estimator);
0661               }
0662             }
0663 
0664             else if (misLayer == (k_LayersAtTIDEnd - 1) ||
0665                      misLayer == k_LayersAtTIDEnd) {  // This is for TID layers 12 and 13
0666               std::vector<ForwardDetLayer const*> negTIDLayers = measTracker.geometricSearchTracker()->negTidLayers();
0667               std::vector<ForwardDetLayer const*> posTIDLayers = measTracker.geometricSearchTracker()->posTidLayers();
0668               const DetLayer* tidLayerneg = negTIDLayers[misLayer - k_LayersAtTOBEnd - 1];
0669               const DetLayer* tidLayerpos = posTIDLayers[misLayer - k_LayersAtTOBEnd - 1];
0670 
0671               if (tTopo->tidSide(iidd) == 1) {
0672                 tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, prop, chi2Estimator);
0673               } else if (tTopo->tidSide(iidd) == 2) {
0674                 tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, prop, chi2Estimator);
0675               }
0676             }
0677 
0678             if (misLayer > k_LayersStart && misLayer < k_LayersAtTOBEnd) {  // Barrel
0679 
0680               std::vector<BarrelDetLayer const*> barrelTIBLayers = measTracker.geometricSearchTracker()->tibLayers();
0681               std::vector<BarrelDetLayer const*> barrelTOBLayers = measTracker.geometricSearchTracker()->tobLayers();
0682 
0683               if (misLayer > k_LayersStart && misLayer <= k_LayersAtTIBEnd) {
0684                 const DetLayer* tibLayer = barrelTIBLayers[misLayer - k_LayersStart - 1];
0685                 tmpTmeas = layerMeasurements.measurements(*tibLayer, tsos, prop, chi2Estimator);
0686               } else if (misLayer > k_LayersAtTIBEnd && misLayer < k_LayersAtTOBEnd) {
0687                 const DetLayer* tobLayer = barrelTOBLayers[misLayer - k_LayersAtTIBEnd - 1];
0688                 tmpTmeas = layerMeasurements.measurements(*tobLayer, tsos, prop, chi2Estimator);
0689               }
0690             }
0691           }
0692           if ((int(TKlayers) > k_LayersStart && int(TKlayers) <= k_LayersAtTIBEnd) && int(prev_TKlayers) == 12) {
0693             const DetLayer* detlayer = itm->layer();
0694             const LayerMeasurements layerMeasurements{measTracker, *measurementTrackerEvent};
0695             const TrajectoryStateOnSurface tsos = itm->updatedState();
0696             std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, prop, chi2Estimator);
0697             std::vector<ForwardDetLayer const*> negTIDLayers = measTracker.geometricSearchTracker()->negTidLayers();
0698             std::vector<ForwardDetLayer const*> posTIDLayers = measTracker.geometricSearchTracker()->posTidLayers();
0699 
0700             const DetLayer* tidLayerneg = negTIDLayers[k_LayersStart];
0701             const DetLayer* tidLayerpos = posTIDLayers[k_LayersStart];
0702             if (tTopo->tidSide(iidd) == 1) {
0703               tmpTmeas = layerMeasurements.measurements(*tidLayerneg, tsos, prop, chi2Estimator);
0704             } else if (tTopo->tidSide(iidd) == 2) {
0705               tmpTmeas = layerMeasurements.measurements(*tidLayerpos, tsos, prop, chi2Estimator);
0706             }
0707           }
0708 
0709           if ((int(TKlayers) > k_LayersAtTIBEnd && int(TKlayers) <= k_LayersAtTOBEnd) && int(prev_TKlayers) == 15) {
0710             const DetLayer* detlayer = itm->layer();
0711             const LayerMeasurements layerMeasurements{measTracker, *measurementTrackerEvent};
0712             const TrajectoryStateOnSurface tsos = itm->updatedState();
0713             std::vector<DetLayer::DetWithState> compatDets = detlayer->compatibleDets(tsos, prop, chi2Estimator);
0714 
0715             std::vector<ForwardDetLayer const*> negTECLayers = measTracker.geometricSearchTracker()->negTecLayers();
0716             std::vector<ForwardDetLayer const*> posTECLayers = measTracker.geometricSearchTracker()->posTecLayers();
0717 
0718             const DetLayer* tecLayerneg = negTECLayers[k_LayersStart];
0719             const DetLayer* tecLayerpos = posTECLayers[k_LayersStart];
0720             if (tTopo->tecSide(iidd) == 1) {
0721               tmpTmeas = layerMeasurements.measurements(*tecLayerneg, tsos, prop, chi2Estimator);
0722             } else if (tTopo->tecSide(iidd) == 2) {
0723               tmpTmeas = layerMeasurements.measurements(*tecLayerpos, tsos, prop, chi2Estimator);
0724             }
0725           }
0726 
0727           if (!tmpTmeas.empty()) {
0728             TrajectoryMeasurement TM_tmp(tmpTmeas.back());
0729             unsigned int iidd_tmp = TM_tmp.recHit()->geographicalId().rawId();
0730             if (iidd_tmp != 0) {
0731               LogDebug("SiStripHitEfficiency:HitEff") << " hit actually being added to TM vector";
0732               if ((!useAllHitsFromTracksWithMissingHits_ || (!useFirstMeas_ && isFirstMeas)))
0733                 TMs.clear();
0734               if (::isDoubleSided(iidd_tmp, tTopo)) {
0735                 TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator, 1));
0736                 TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator, 2));
0737               } else
0738                 TMs.push_back(TrajectoryAtInvalidHit(TM_tmp, tTopo, tkgeom, propagator));
0739               missingHitAdded = true;
0740               hitRecoveryCounters[misLayer] += 1;
0741             }
0742           }
0743         }
0744 
0745         prev_TKlayers = TKlayers;
0746         if (!useFirstMeas_ && isFirstMeas && !missingHitAdded)
0747           continue;
0748         if (!useLastMeas_ && isLastMeas)
0749           continue;
0750         bool hitsWithBias = false;
0751         for (auto ilayer : missedLayers) {
0752           if (ilayer < TKlayers)
0753             hitsWithBias = true;
0754         }
0755         if (hasMissingHits && theInHit->getType() != TrackingRecHit::Type::missing && !missingHitAdded &&
0756             hitsWithBias && !useAllHitsFromTracksWithMissingHits_) {
0757           continue;
0758         }
0759 
0760         //////////////////////////////////////////////
0761         //Now check for tracks at TOB6 and TEC9
0762 
0763         // to make sure we only propagate on the last TOB5 hit check the next entry isn't also in TOB5
0764         // to avoid bias, make sure the TOB5 hit is valid (an invalid hit on TOB5 could only exist with a valid hit on TOB6)
0765         const auto nextId = (itm + 1 != TMeas.end()) ? (itm + 1)->recHit()->geographicalId() : DetId{};  // null if last
0766 
0767         if (TKlayers == 9 && theInHit->isValid() && !((!nextId.null()) && (::checkLayer(nextId.rawId(), tTopo) == 9))) {
0768           //      if ( TKlayers==9 && itm==TMeas.rbegin()) {
0769           //      if ( TKlayers==9 && (itm==TMeas.back()) ) {     // to check for only the last entry in the trajectory for propagation
0770           const DetLayer* tob6 = measTracker.geometricSearchTracker()->tobLayers().back();
0771           const LayerMeasurements theLayerMeasurements{measTracker, *measurementTrackerEvent};
0772           const TrajectoryStateOnSurface tsosTOB5 = itm->updatedState();
0773           const auto tmp = theLayerMeasurements.measurements(*tob6, tsosTOB5, prop, chi2Estimator);
0774 
0775           if (!tmp.empty()) {
0776             LogDebug("SiStripHitEfficiencyWorker") << "size of TM from propagation = " << tmp.size();
0777 
0778             // take the last of the TMs, which is always an invalid hit
0779             // if no detId is available, ie detId==0, then no compatible layer was crossed
0780             // otherwise, use that TM for the efficiency measurement
0781             const auto& tob6TM = tmp.back();
0782             const auto& tob6Hit = tob6TM.recHit();
0783             if (tob6Hit->geographicalId().rawId() != 0) {
0784               LogDebug("SiStripHitEfficiencyWorker") << "tob6 hit actually being added to TM vector";
0785               TMs.emplace_back(tob6TM, tTopo, tkgeom, propagator);
0786             }
0787           }
0788         }
0789 
0790         // same for TEC8
0791         if (TKlayers == 21 && theInHit->isValid() &&
0792             !((!nextId.null()) && (::checkLayer(nextId.rawId(), tTopo) == 21))) {
0793           const DetLayer* tec9pos = measTracker.geometricSearchTracker()->posTecLayers().back();
0794           const DetLayer* tec9neg = measTracker.geometricSearchTracker()->negTecLayers().back();
0795 
0796           const LayerMeasurements theLayerMeasurements{measTracker, *measurementTrackerEvent};
0797           const TrajectoryStateOnSurface tsosTEC9 = itm->updatedState();
0798 
0799           // check if track on positive or negative z
0800           if (!(iidd == SiStripSubdetector::TEC))
0801             LogDebug("SiStripHitEfficiencyWorker") << "there is a problem with TEC 9 extrapolation";
0802 
0803           //LogDebug("SiStripHitEfficiencyWorker") << " tec9 id = " << iidd << " and side = " << tTopo->tecSide(iidd) ;
0804           std::vector<TrajectoryMeasurement> tmp;
0805           if (tTopo->tecSide(iidd) == 1) {
0806             tmp = theLayerMeasurements.measurements(*tec9neg, tsosTEC9, prop, chi2Estimator);
0807             //LogDebug("SiStripHitEfficiencyWorker") << "on negative side" ;
0808           }
0809           if (tTopo->tecSide(iidd) == 2) {
0810             tmp = theLayerMeasurements.measurements(*tec9pos, tsosTEC9, prop, chi2Estimator);
0811             //LogDebug("SiStripHitEfficiencyWorker") << "on positive side" ;
0812           }
0813 
0814           if (!tmp.empty()) {
0815             // take the last of the TMs, which is always an invalid hit
0816             // if no detId is available, ie detId==0, then no compatible layer was crossed
0817             // otherwise, use that TM for the efficiency measurement
0818             const auto& tec9TM = tmp.back();
0819             const auto& tec9Hit = tec9TM.recHit();
0820 
0821             const unsigned int tec9id = tec9Hit->geographicalId().rawId();
0822             LogDebug("SiStripHitEfficiencyWorker")
0823                 << "tec9id = " << tec9id << " is Double sided = " << ::isDoubleSided(tec9id, tTopo)
0824                 << "  and 0x3 = " << (tec9id & 0x3);
0825 
0826             if (tec9Hit->geographicalId().rawId() != 0) {
0827               LogDebug("SiStripHitEfficiencyWorker") << "tec9 hit actually being added to TM vector";
0828               // in tec the hit can be single or doubled sided. whenever the invalid hit at the end of vector of TMs is
0829               // double sided it is always on the matched surface, so we need to split it into the true sensor surfaces
0830               if (::isDoubleSided(tec9id, tTopo)) {
0831                 TMs.emplace_back(tec9TM, tTopo, tkgeom, propagator, 1);
0832                 TMs.emplace_back(tec9TM, tTopo, tkgeom, propagator, 2);
0833               } else
0834                 TMs.emplace_back(tec9TM, tTopo, tkgeom, propagator);
0835             }
0836           }  //else LogDebug("SiStripHitEfficiencyWorker") << "tec9 tmp empty" ;
0837         }
0838 
0839         for (const auto& tm : TMs) {
0840           fillForTraj(tm,
0841                       tTopo,
0842                       tkgeom,
0843                       stripcpe,
0844                       stripQuality,
0845                       fedErrorIds,
0846                       commonModeDigis,
0847                       *theClusters,
0848                       e.bunchCrossing(),
0849                       instLumi,
0850                       PU,
0851                       highPurity);
0852         }
0853         LogDebug("SiStripHitEfficiencyWorker") << "After looping over TrajAtValidHit list";
0854       }
0855       LogDebug("SiStripHitEfficiencyWorker") << "end TMeasurement loop";
0856     }
0857     LogDebug("SiStripHitEfficiencyWorker") << "end of trajectories loop";
0858   }
0859 }
0860 
0861 void SiStripHitEfficiencyWorker::fillForTraj(const TrajectoryAtInvalidHit& tm,
0862                                              const TrackerTopology* tTopo,
0863                                              const TrackerGeometry* tkgeom,
0864                                              const StripClusterParameterEstimator& stripCPE,
0865                                              const SiStripQuality& stripQuality,
0866                                              const DetIdVector& fedErrorIds,
0867                                              const edm::Handle<edm::DetSetVector<SiStripRawDigi>>& commonModeDigis,
0868                                              const edmNew::DetSetVector<SiStripCluster>& theClusters,
0869                                              int bunchCrossing,
0870                                              float instLumi,
0871                                              float PU,
0872                                              bool highPurity) {
0873   // --> Get trajectory from combinatedStat& e
0874   const auto iidd = tm.monodet_id();
0875   LogDebug("SiStripHitEfficiencyWorker") << "setting iidd = " << iidd << " before checking efficiency and ";
0876 
0877   const auto xloc = tm.localX();
0878   const auto yloc = tm.localY();
0879 
0880   const auto xErr = tm.localErrorX();
0881   const auto yErr = tm.localErrorY();
0882 
0883   int TrajStrip = -1;
0884 
0885   // reget layer from iidd here, to account for TOB 6 and TEC 9 TKlayers being off
0886   const auto TKlayers = ::checkLayer(iidd, tTopo);
0887 
0888   const bool withinAcceptance =
0889       tm.withinAcceptance() && (!::isInBondingExclusionZone(iidd, TKlayers, yloc, yErr, tTopo));
0890 
0891   if (  // (TKlayers > 0) && // FIXME confirm this
0892       ((layers_ == TKlayers) ||
0893        (layers_ == bounds::k_LayersStart))) {  // Look at the layer not used to reconstruct the track
0894     LogDebug("SiStripHitEfficiencyWorker") << "Looking at layer under study";
0895     unsigned int ModIsBad = 2;
0896     unsigned int SiStripQualBad = 0;
0897     float commonMode = -100;
0898 
0899     // RPhi RecHit Efficiency
0900 
0901     if (!theClusters.empty()) {
0902       LogDebug("SiStripHitEfficiencyWorker") << "Checking clusters with size = " << theClusters.size();
0903       std::vector<::ClusterInfo> VCluster_info;  //fill with X residual, X residual pull, local X
0904       const auto idsv = theClusters.find(iidd);
0905       if (idsv != theClusters.end()) {
0906         //if (DEBUG_)      LogDebug("SiStripHitEfficiencyWorker") << "the ID from the dsv = " << dsv.id() ;
0907         LogDebug("SiStripHitEfficiencyWorker")
0908             << "found  (ClusterId == iidd) with ClusterId = " << idsv->id() << " and iidd = " << iidd;
0909         const auto stripdet = dynamic_cast<const StripGeomDetUnit*>(tkgeom->idToDetUnit(DetId(iidd)));
0910         const StripTopology& Topo = stripdet->specificTopology();
0911 
0912         float hbedge = 0.0;
0913         float htedge = 0.0;
0914         float hapoth = 0.0;
0915         float uylfac = 0.0;
0916         float uxlden = 0.0;
0917         if (TKlayers > bounds::k_LayersAtTOBEnd) {
0918           const BoundPlane& plane = stripdet->surface();
0919           const TrapezoidalPlaneBounds* trapezoidalBounds(
0920               dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
0921           std::array<const float, 4> const& parameterTrap = (*trapezoidalBounds).parameters();  // el bueno aqui
0922           hbedge = parameterTrap[0];
0923           htedge = parameterTrap[1];
0924           hapoth = parameterTrap[3];
0925           uylfac = (htedge - hbedge) / (htedge + hbedge) / hapoth;
0926           uxlden = 1 + yloc * uylfac;
0927         }
0928 
0929         // Need to know position of trajectory in strip number for selecting the right APV later
0930         if (TrajStrip == -1) {
0931           int nstrips = Topo.nstrips();
0932           float pitch = stripdet->surface().bounds().width() / nstrips;
0933           TrajStrip = xloc / pitch + nstrips / 2.0;
0934           // Need additionnal corrections for endcap
0935           if (TKlayers > bounds::k_LayersAtTOBEnd) {
0936             const float TrajLocXMid = xloc / (1 + (htedge - hbedge) * yloc / (htedge + hbedge) /
0937                                                       hapoth);  // radialy extrapolated x loc position at middle
0938             TrajStrip = TrajLocXMid / pitch + nstrips / 2.0;
0939           }
0940           //LogDebug("SiStripHitEfficiency")<<" Layer "<<TKlayers<<" TrajStrip: "<<nstrips<<" "<<pitch<<" "<<TrajStrip;;
0941         }
0942 
0943         for (const auto& clus : *idsv) {
0944           StripClusterParameterEstimator::LocalValues parameters = stripCPE.localParameters(clus, *stripdet);
0945           float res = (parameters.first.x() - xloc);
0946           float sigma = ::checkConsistency(parameters, xloc, xErr);
0947           // The consistency is probably more accurately measured with the Chi2MeasurementEstimator. To use it
0948           // you need a TransientTrackingRecHit instead of the cluster
0949           //theEstimator=       new Chi2MeasurementEstimator(30);
0950           //const Chi2MeasurementEstimator *theEstimator(100);
0951           //theEstimator->estimate(tm.tsos(), TransientTrackingRecHit);
0952 
0953           if (TKlayers > bounds::k_LayersAtTOBEnd) {
0954             res = parameters.first.x() - xloc / uxlden;  // radialy extrapolated x loc position at middle
0955             sigma = abs(res) / sqrt(parameters.second.xx() + xErr * xErr / uxlden / uxlden +
0956                                     yErr * yErr * xloc * xloc * uylfac * uylfac / uxlden / uxlden / uxlden / uxlden);
0957           }
0958 
0959           VCluster_info.emplace_back(res, sigma, parameters.first.x());
0960 
0961           LogDebug("SiStripHitEfficiencyWorker") << "Have ID match. residual = " << res << "  res sigma = " << sigma;
0962           //LogDebug("SiStripHitEfficiencyWorker")
0963           //    << "trajectory measurement compatability estimate = " << (*itm).estimate() ;
0964           LogDebug("SiStripHitEfficiencyWorker")
0965               << "hit position = " << parameters.first.x() << "  hit error = " << sqrt(parameters.second.xx())
0966               << "  trajectory position = " << xloc << "  traj error = " << xErr;
0967         }
0968       }
0969       ::ClusterInfo finalCluster{1000.0, 1000.0, 0.0};
0970       if (!VCluster_info.empty()) {
0971         LogDebug("SiStripHitEfficiencyWorker") << "found clusters > 0";
0972         if (VCluster_info.size() > 1) {
0973           //get the smallest one
0974           for (const auto& res : VCluster_info) {
0975             if (std::abs(res.xResidualPull) < std::abs(finalCluster.xResidualPull)) {
0976               finalCluster = res;
0977             }
0978             LogDebug("SiStripHitEfficiencyWorker")
0979                 << "iresidual = " << res.xResidual << "  isigma = " << res.xResidualPull
0980                 << "  and FinalRes = " << finalCluster.xResidual;
0981           }
0982         } else {
0983           finalCluster = VCluster_info[0];
0984         }
0985         VCluster_info.clear();
0986       }
0987 
0988       LogDebug("SiStripHitEfficiencyWorker") << "Final residual in X = " << finalCluster.xResidual << "+-"
0989                                              << (finalCluster.xResidual / finalCluster.xResidualPull);
0990       LogDebug("SiStripHitEfficiencyWorker")
0991           << "Checking location of trajectory: abs(yloc) = " << abs(yloc) << "  abs(xloc) = " << abs(xloc);
0992 
0993       //
0994       // fill ntuple varibles
0995 
0996       //if ( stripQuality->IsModuleBad(iidd) )
0997       if (stripQuality.getBadApvs(iidd) != 0) {
0998         SiStripQualBad = 1;
0999         LogDebug("SiStripHitEfficiencyWorker") << "strip is bad from SiStripQuality";
1000       } else {
1001         SiStripQualBad = 0;
1002         LogDebug("SiStripHitEfficiencyWorker") << "strip is good from SiStripQuality";
1003       }
1004 
1005       //check for FED-detected errors and include those in SiStripQualBad
1006       for (unsigned int ii = 0; ii < fedErrorIds.size(); ii++) {
1007         if (iidd == fedErrorIds[ii].rawId())
1008           SiStripQualBad = 1;
1009       }
1010 
1011       // CM of APV crossed by traj
1012       if (addCommonMode_)
1013         if (commonModeDigis.isValid() && TrajStrip >= 0 && TrajStrip <= 768) {
1014           const auto digiframe = commonModeDigis->find(iidd);
1015           if (digiframe != commonModeDigis->end())
1016             if ((unsigned)TrajStrip / sistrip::STRIPS_PER_APV < digiframe->data.size())
1017               commonMode = digiframe->data.at(TrajStrip / sistrip::STRIPS_PER_APV).adc();
1018         }
1019 
1020       LogDebug("SiStripHitEfficiencyWorker") << "before check good";
1021 
1022       if (finalCluster.xResidualPull < 999.0) {  //could make requirement on track/hit consistency, but for
1023         //now take anything with a hit on the module
1024         LogDebug("SiStripHitEfficiencyWorker")
1025             << "hit being counted as good " << finalCluster.xResidual << " FinalRecHit " << iidd << "   TKlayers  "
1026             << TKlayers << " xloc " << xloc << " yloc  " << yloc << " module " << iidd
1027             << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/" << ((iidd & 0x3) == 1) << "/"
1028             << ((iidd & 0x3) == 2);
1029         ModIsBad = 0;
1030       } else {
1031         LogDebug("SiStripHitEfficiencyWorker")
1032             << "hit being counted as bad   ######### Invalid RPhi FinalResX " << finalCluster.xResidual
1033             << " FinalRecHit " << iidd << "   TKlayers  " << TKlayers << " xloc " << xloc << " yloc  " << yloc
1034             << " module " << iidd << "   matched/stereo/rphi = " << ((iidd & 0x3) == 0) << "/" << ((iidd & 0x3) == 1)
1035             << "/" << ((iidd & 0x3) == 2);
1036         ModIsBad = 1;
1037         LogDebug("SiStripHitEfficiencyWorker")
1038             << " RPhi Error " << sqrt(xErr * xErr + yErr * yErr) << " ErrorX " << xErr << " yErr " << yErr;
1039       }
1040 
1041       LogDebug("SiStripHitEfficiencyWorker")
1042           << "To avoid them staying unused: ModIsBad=" << ModIsBad << ", SiStripQualBad=" << SiStripQualBad
1043           << ", commonMode=" << commonMode << ", highPurity=" << highPurity
1044           << ", withinAcceptance=" << withinAcceptance;
1045 
1046       unsigned int layer = TKlayers;
1047       if (showRings_ && layer > bounds::k_LayersAtTOBEnd) {  // use rings instead of wheels
1048         if (layer <= bounds::k_LayersAtTIDEnd) {             // TID
1049           layer = bounds::k_LayersAtTOBEnd +
1050                   tTopo->tidRing(iidd);  // ((iidd >> 9) & 0x3);  // 3 disks and also 3 rings -> use the same container
1051         } else {                         // TEC
1052           layer = bounds::k_LayersAtTIDEnd + tTopo->tecRing(iidd);  // ((iidd >> 5) & 0x7);
1053         }
1054       }
1055       unsigned int layerWithSide = layer;
1056       if (layer > bounds::k_LayersAtTOBEnd && layer <= bounds::k_LayersAtTIDEnd) {
1057         const auto side = tTopo->tidSide(iidd);  //(iidd >> 13) & 0x3;  // TID
1058         if (side == 2)
1059           layerWithSide = layer + 3;
1060       } else if (layer > bounds::k_LayersAtTIDEnd) {
1061         const auto side = tTopo->tecSide(iidd);  // (iidd >> 18) & 0x3;  // TEC
1062         if (side == 1) {
1063           layerWithSide = layer + 3;
1064         } else if (side == 2) {
1065           layerWithSide = layer + 3 + (showRings_ ? 7 : 9);
1066         }
1067       }
1068 
1069       if ((bunchX_ > 0 && bunchX_ != bunchCrossing) || (!withinAcceptance) ||
1070           (useOnlyHighPurityTracks_ && !highPurity) ||
1071           (!showTOB6TEC9_ && (TKlayers == bounds::k_LayersAtTOBEnd || TKlayers == bounds::k_LayersAtTECEnd)) ||
1072           (badModules_.end() != badModules_.find(iidd)))
1073         return;
1074 
1075       const bool badquality = (SiStripQualBad == 1);
1076 
1077       //Now that we have a good event, we need to look at if we expected it or not, and the location
1078       //if we didn't
1079       //Fill the missing hit information first
1080       bool badflag = false;  // true for hits that are expected but not found
1081       if (resXSig_ < 0) {
1082         if (ModIsBad == 1)
1083           badflag = true;  // isBad set to false in the tree when resxsig<999.0
1084       } else {
1085         if (ModIsBad == 1 || finalCluster.xResidualPull > resXSig_)
1086           badflag = true;
1087       }
1088 
1089       // Conversion of positions in strip unit
1090       int nstrips = -9;
1091       float Pitch = -9.0;
1092       const StripGeomDetUnit* stripdet = nullptr;
1093       if (finalCluster.xResidualPull ==
1094           1000.0) {      // special treatment, no GeomDetUnit associated in some cases when no cluster found
1095         Pitch = 0.0205;  // maximum
1096         nstrips = 768;   // maximum
1097       } else {
1098         stripdet = dynamic_cast<const StripGeomDetUnit*>(tkgeom->idToDetUnit(iidd));
1099         const StripTopology& Topo = stripdet->specificTopology();
1100         nstrips = Topo.nstrips();
1101         Pitch = stripdet->surface().bounds().width() / Topo.nstrips();
1102       }
1103       double stripTrajMid = xloc / Pitch + nstrips / 2.0;
1104       double stripCluster = finalCluster.xLocal / Pitch + nstrips / 2.0;
1105       // For trapezoidal modules: extrapolation of x trajectory position to the y middle of the module
1106       //  for correct comparison with cluster position
1107       if (stripdet && layer > bounds::k_LayersAtTOBEnd) {
1108         const auto& trapezoidalBounds = dynamic_cast<const TrapezoidalPlaneBounds&>(stripdet->surface().bounds());
1109         std::array<const float, 4> const& parameters = trapezoidalBounds.parameters();
1110         const float hbedge = parameters[0];
1111         const float htedge = parameters[1];
1112         const float hapoth = parameters[3];
1113         const float TrajLocXMid = xloc / (1 + (htedge - hbedge) * yloc / (htedge + hbedge) /
1114                                                   hapoth);  // radialy extrapolated x loc position at middle
1115         stripTrajMid = TrajLocXMid / Pitch + nstrips / 2.0;
1116       }
1117 
1118       if ((!badquality) && (layer < h_resolution.size())) {
1119         LogDebug("SiStripHitEfficiencyWorker")
1120             << "layer " << layer << " vector index " << layer - 1 << " before filling h_resolution" << std::endl;
1121         h_resolution[layer - 1]->Fill(finalCluster.xResidualPull != 1000.0 ? stripTrajMid - stripCluster : 1000);
1122       }
1123 
1124       // New matching methods
1125       if (clusterMatchingMethod_ >= 1) {
1126         badflag = false;
1127         if (finalCluster.xResidualPull == 1000.0) {
1128           LogDebug("SiStripHitEfficiencyWorker") << "Marking bad for resxsig=1000";
1129           badflag = true;
1130         } else {
1131           if (clusterMatchingMethod_ == 2 || clusterMatchingMethod_ == 4) {
1132             // check the distance between cluster and trajectory position
1133             if (std::abs(stripCluster - stripTrajMid) > clusterTracjDist_) {
1134               LogDebug("SiStripHitEfficiencyWorker") << "Marking bad for cluster-to-traj distance";
1135               badflag = true;
1136             }
1137           }
1138           if (clusterMatchingMethod_ == 3 || clusterMatchingMethod_ == 4) {
1139             // cluster and traj have to be in the same APV (don't take edges into accounts)
1140             const int tapv = (int)stripTrajMid / sistrip::STRIPS_PER_APV;
1141             const int capv = (int)stripCluster / sistrip::STRIPS_PER_APV;
1142             float stripInAPV = stripTrajMid - tapv * sistrip::STRIPS_PER_APV;
1143             if (stripInAPV < stripsApvEdge_ || stripInAPV > sistrip::STRIPS_PER_APV - stripsApvEdge_) {
1144               LogDebug("SiStripHitEfficiencyWorker") << "Too close to the edge: " << stripInAPV;
1145               return;
1146             }
1147             if (tapv != capv) {
1148               LogDebug("SiStripHitEfficiencyWorker") << "Marking bad for tapv!=capv";
1149               badflag = true;
1150             }
1151           }
1152         }
1153       }
1154       if (!badquality) {
1155         LogDebug("SiStripHitEfficiencyWorker")
1156             << "Filling measurement for " << iidd << " in layer " << layer << " histograms with bx=" << bunchCrossing
1157             << ", lumi=" << instLumi << ", PU=" << PU << "; bad flag=" << badflag;
1158 
1159         // hot/cold maps of hits that are expected but not found
1160         if (badflag) {
1161           if (layer > bounds::k_LayersStart && layer <= bounds::k_LayersAtTIBEnd) {
1162             //We are in the TIB
1163             float phi = ::calcPhi(tm.globalX(), tm.globalY());
1164             h_hotcold[layer - 1]->Fill(360. - phi, tm.globalZ(), 1.);
1165           } else if (layer > bounds::k_LayersAtTIBEnd && layer <= bounds::k_LayersAtTOBEnd) {
1166             //We are in the TOB
1167             float phi = ::calcPhi(tm.globalX(), tm.globalY());
1168             h_hotcold[layer - 1]->Fill(360. - phi, tm.globalZ(), 1.);
1169           } else if (layer > bounds::k_LayersAtTOBEnd && layer <= bounds::k_LayersAtTIDEnd) {
1170             //We are in the TID
1171             //There are 2 different maps here
1172             int side = tTopo->tidSide(iidd);
1173             if (side == 1)
1174               h_hotcold[(layer - 1) + (layer - 11)]->Fill(-tm.globalY(), tm.globalX(), 1.);
1175             else if (side == 2)
1176               h_hotcold[(layer - 1) + (layer - 10)]->Fill(-tm.globalY(), tm.globalX(), 1.);
1177           } else if (layer > bounds::k_LayersAtTIDEnd) {
1178             //We are in the TEC
1179             //There are 2 different maps here
1180             int side = tTopo->tecSide(iidd);
1181             if (side == 1)
1182               h_hotcold[(layer + 2) + (layer - 14)]->Fill(-tm.globalY(), tm.globalX(), 1.);
1183             else if (side == 2)
1184               h_hotcold[(layer + 2) + (layer - 13)]->Fill(-tm.globalY(), tm.globalX(), 1.);
1185           }
1186         }
1187 
1188         LogDebug("SiStripHitEfficiencyWorker")
1189             << "layer " << layer << " vector index " << layer - 1 << " before filling h_layer_vsSmthg" << std::endl;
1190         h_layer_vsBx[layer - 1].fill(bunchCrossing, !badflag);
1191         if (addLumi_) {
1192           h_layer_vsLumi[layer - 1].fill(instLumi, !badflag);
1193           h_layer_vsPU[layer - 1].fill(PU, !badflag);
1194         }
1195         if (addCommonMode_) {
1196           h_layer_vsCM[layer - 1].fill(commonMode, !badflag);
1197         }
1198         h_goodLayer.fill(layerWithSide, !badflag);
1199       }
1200       // efficiency without bad modules excluded
1201       h_allLayer.fill(layerWithSide, !badflag);
1202 
1203       // efficiency without bad modules excluded
1204       if (TKlayers) {
1205         h_module.fill(iidd, !badflag);
1206         assert(h_module.check(iidd));
1207       }
1208 
1209       /* Used in SiStripHitEffFromCalibTree:
1210        * run              -> "run"              -> run              // e.id().run()
1211        * event            -> "event"            -> evt              // e.id().event()
1212        * ModIsBad         -> "ModIsBad"         -> isBad
1213        * SiStripQualBad   -> "SiStripQualBad""  -> quality
1214        * Id               -> "Id"               -> id               // iidd
1215        * withinAcceptance -> "withinAcceptance" -> accept
1216        * whatlayer        -> "layer"            -> layer_wheel      // Tklayers
1217        * highPurity       -> "highPurity"       -> highPurity
1218        * TrajGlbX         -> "TrajGlbX"         -> x                // tm.globalX()
1219        * TrajGlbY         -> "TrajGlbY"         -> y                // tm.globalY()
1220        * TrajGlbZ         -> "TrajGlbZ"         -> z                // tm.globalZ()
1221        * ResXSig          -> "ResXSig"          -> resxsig          // finalCluster.xResidualPull;
1222        * TrajLocX         -> "TrajLocX"         -> TrajLocX         // xloc
1223        * TrajLocY         -> "TrajLocY"         -> TrajLocY         // yloc
1224        * ClusterLocX      -> "ClusterLocX"      -> ClusterLocX      // finalCluster.xLocal
1225        * bunchx           -> "bunchx"           -> bx               // e.bunchCrossing()
1226        * instLumi         -> "instLumi"         -> instLumi         ## if addLumi_
1227        * PU               -> "PU"               -> PU               ## if addLumi_
1228        * commonMode       -> "commonMode"       -> CM               ## if addCommonMode_ / _useCM
1229        */
1230       LogDebug("SiStripHitEfficiencyWorker") << "after good location check";
1231     }
1232     LogDebug("SiStripHitEfficiencyWorker") << "after list of clusters";
1233   }
1234   LogDebug("SiStripHitEfficiencyWorker") << "After layers=TKLayers if with TKlayers=" << TKlayers
1235                                          << ", layers=" << layers_;
1236 }
1237 
1238 void SiStripHitEfficiencyWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1239   edm::ParameterSetDescription desc;
1240   desc.add<std::string>("dqmDir", "AlCaReco/SiStripHitEfficiency");
1241   desc.add<bool>("UseOnlyHighPurityTracks", true);
1242   desc.add<bool>("cutOnTracks", false);
1243   desc.add<bool>("doMissingHitsRecovery", false);
1244   desc.add<bool>("useAllHitsFromTracksWithMissingHits", false);
1245   desc.add<bool>("useFirstMeas", false);
1246   desc.add<bool>("useLastMeas", false);
1247   desc.add<double>("ClusterTrajDist", 64.0);
1248   desc.add<double>("ResXSig", -1);
1249   desc.add<double>("StripsApvEdge", 10.0);
1250   desc.add<edm::InputTag>("combinatorialTracks", edm::InputTag{"generalTracks"});
1251   desc.add<edm::InputTag>("commonMode", edm::InputTag{"siStripDigis", "CommonMode"});
1252   desc.add<edm::InputTag>("lumiScalers", edm::InputTag{"scalersRawToDigi"});
1253   desc.add<edm::InputTag>("metadata", edm::InputTag{"onlineMetaDataDigis"});
1254   desc.add<edm::InputTag>("siStripClusters", edm::InputTag{"siStripClusters"});
1255   desc.add<edm::InputTag>("siStripDigis", edm::InputTag{"siStripDigis"});
1256   desc.add<edm::InputTag>("trackerEvent", edm::InputTag{"MeasurementTrackerEvent"});
1257   desc.add<edm::InputTag>("trajectories", edm::InputTag{"generalTracks"});
1258   desc.add<int>("ClusterMatchingMethod", 0);
1259   desc.add<int>("Layer", 0);
1260   desc.add<unsigned int>("trackMultiplicity", 100);
1261   desc.addUntracked<bool>("Debug", false);
1262   desc.addUntracked<bool>("ShowRings", false);
1263   desc.addUntracked<bool>("ShowTOB6TEC9", false);
1264   desc.addUntracked<bool>("addCommonMode", false);
1265   desc.addUntracked<bool>("addLumi", true);
1266   desc.addUntracked<int>("BunchCrossing", 0);
1267   desc.addUntracked<std::string>("BadModulesFile", "");
1268   descriptions.addWithDefaultLabel(desc);
1269 }
1270 
1271 #include "FWCore/Framework/interface/MakerMacros.h"
1272 DEFINE_FWK_MODULE(SiStripHitEfficiencyWorker);