Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-02-20 03:45:59

0001 // -*- C++ -*-
0002 //bookLayer
0003 // Package:    Phase2OTValidateCluster
0004 // Class:      Phase2OTValidateCluster
0005 //
0006 /**\class Phase2OTValidateCluster Phase2OTValidateCluster.cc 
0007 
0008  Description: Validation plots tracker clusters. 
0009 
0010 */
0011 //
0012 // Author: Gabriel Ramirez, Suvankar Roy Chowdhury
0013 // Date: May 23, 2020
0014 //
0015 #include <memory>
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/Framework/interface/Event.h"
0019 #include "FWCore/Framework/interface/ESHandle.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 #include "FWCore/Framework/interface/ESWatcher.h"
0022 #include "FWCore/Framework/interface/Frameworkfwd.h"
0023 #include "FWCore/ServiceRegistry/interface/Service.h"
0024 #include "FWCore/Utilities/interface/InputTag.h"
0025 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0026 
0027 #include "DataFormats/Common/interface/Handle.h"
0028 #include "DataFormats/DetId/interface/DetId.h"
0029 #include "DataFormats/Phase2TrackerCluster/interface/Phase2TrackerCluster1D.h"
0030 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0031 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0032 #include "DataFormats/Common/interface/DetSetVector.h"
0033 
0034 #include "SimDataFormats/TrackerDigiSimLink/interface/PixelDigiSimLink.h"
0035 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0036 #include "SimDataFormats/Vertex/interface/SimVertexContainer.h"
0037 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0038 
0039 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0040 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0041 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0042 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0043 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0044 #include "DataFormats/Common/interface/DetSetVector.h"
0045 // DQM Histograming
0046 #include "DQMServices/Core/interface/MonitorElement.h"
0047 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0048 #include "DQMServices/Core/interface/DQMStore.h"
0049 
0050 class Phase2OTValidateCluster : public DQMEDAnalyzer {
0051 public:
0052   typedef std::map<unsigned int, std::vector<PSimHit>> SimHitsMap;
0053   typedef std::map<unsigned int, SimTrack> SimTracksMap;
0054 
0055   explicit Phase2OTValidateCluster(const edm::ParameterSet&);
0056   ~Phase2OTValidateCluster() override;
0057   void bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) override;
0058   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0059   void dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) override;
0060   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0061 
0062 private:
0063   struct ClusterMEs {
0064     MonitorElement* deltaX_S = nullptr;
0065     MonitorElement* deltaX_P = nullptr;
0066     MonitorElement* deltaY_S = nullptr;
0067     MonitorElement* deltaY_P = nullptr;
0068     MonitorElement* deltaX_P_primary = nullptr;
0069     MonitorElement* deltaY_P_primary = nullptr;
0070     MonitorElement* deltaX_S_primary = nullptr;
0071     MonitorElement* deltaY_S_primary = nullptr;
0072   };
0073 
0074   void fillOTHistos(const edm::Event& iEvent,
0075                     const std::vector<edm::Handle<edm::PSimHitContainer>>& simHits,
0076                     const std::map<unsigned int, SimTrack>& simTracks);
0077   void bookLayerHistos(DQMStore::IBooker& ibooker, uint32_t det_it, const std::string& subdir);
0078   std::vector<unsigned int> getSimTrackId(const edm::Handle<edm::DetSetVector<PixelDigiSimLink>>& pixelSimLinks,
0079                                           const DetId& detId,
0080                                           unsigned int channel);
0081 
0082   std::map<std::string, ClusterMEs> layerMEs_;
0083 
0084   edm::ParameterSet config_;
0085   double simtrackminpt_;
0086   std::vector<edm::EDGetTokenT<edm::PSimHitContainer>> simHitTokens_;
0087   edm::EDGetTokenT<edm::DetSetVector<PixelDigiSimLink>> simOTLinksToken_;
0088   edm::EDGetTokenT<edm::DetSetVector<PixelDigiSimLink>> simITLinksToken_;
0089   edm::EDGetTokenT<edm::SimTrackContainer> simTracksToken_;
0090   edm::EDGetTokenT<Phase2TrackerCluster1DCollectionNew> clustersToken_;
0091   std::vector<edm::InputTag> pSimHitSrc_;
0092   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0093   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoToken_;
0094   const TrackerGeometry* tkGeom_ = nullptr;
0095   const TrackerTopology* tTopo_ = nullptr;
0096 };
0097 #include "Validation/SiTrackerPhase2V/interface/TrackerPhase2ValidationUtil.h"
0098 #include "DQM/SiTrackerPhase2/interface/TrackerPhase2DQMUtil.h"
0099 //
0100 // constructors
0101 //
0102 Phase2OTValidateCluster::Phase2OTValidateCluster(const edm::ParameterSet& iConfig)
0103     : config_(iConfig),
0104       simtrackminpt_(config_.getParameter<double>("SimTrackMinPt")),
0105       simOTLinksToken_(consumes<edm::DetSetVector<PixelDigiSimLink>>(
0106           config_.getParameter<edm::InputTag>("OuterTrackerDigiSimLinkSource"))),
0107       simTracksToken_(consumes<edm::SimTrackContainer>(config_.getParameter<edm::InputTag>("simtracks"))),
0108       clustersToken_(
0109           consumes<Phase2TrackerCluster1DCollectionNew>(config_.getParameter<edm::InputTag>("ClusterSource"))),
0110       pSimHitSrc_(config_.getParameter<std::vector<edm::InputTag>>("PSimHitSource")),
0111       geomToken_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord, edm::Transition::BeginRun>()),
0112       topoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd, edm::Transition::BeginRun>()) {
0113   edm::LogInfo("Phase2OTValidateCluster") << ">>> Construct Phase2OTValidateCluster ";
0114   for (const auto& itag : pSimHitSrc_)
0115     simHitTokens_.push_back(consumes<edm::PSimHitContainer>(itag));
0116 }
0117 
0118 Phase2OTValidateCluster::~Phase2OTValidateCluster() {
0119   // do anything here that needs to be done at desctruction time
0120   // (e.g. close files, deallocate resources etc.)
0121   edm::LogInfo("Phase2OTValidateCluster") << ">>> Destroy Phase2OTValidateCluster ";
0122 }
0123 //
0124 // -- DQM Begin Run
0125 //
0126 void Phase2OTValidateCluster::dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0127   tkGeom_ = &iSetup.getData(geomToken_);
0128   tTopo_ = &iSetup.getData(topoToken_);
0129 }
0130 
0131 // -- Analyze
0132 //
0133 void Phase2OTValidateCluster::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0134   // Getting simHits
0135   std::vector<edm::Handle<edm::PSimHitContainer>> simHits;
0136   for (const auto& itoken : simHitTokens_) {
0137     const auto& simHitHandle = iEvent.getHandle(itoken);
0138     if (!simHitHandle.isValid())
0139       continue;
0140     simHits.emplace_back(simHitHandle);
0141   }
0142   // Get the SimTracks
0143   const auto& simTracksRaw = iEvent.getHandle(simTracksToken_);
0144   // Rearrange the simTracks for ease of use <simTrackID, simTrack>
0145   SimTracksMap simTracks;
0146   for (edm::SimTrackContainer::const_iterator simTrackIt(simTracksRaw->begin()); simTrackIt != simTracksRaw->end();
0147        ++simTrackIt) {
0148     if (simTrackIt->momentum().pt() > simtrackminpt_) {
0149       simTracks.emplace(simTrackIt->trackId(), *simTrackIt);
0150     }
0151   }
0152   fillOTHistos(iEvent, simHits, simTracks);
0153 }
0154 
0155 void Phase2OTValidateCluster::fillOTHistos(const edm::Event& iEvent,
0156                                            const std::vector<edm::Handle<edm::PSimHitContainer>>& simHits,
0157                                            const std::map<unsigned int, SimTrack>& simTracks) {
0158   // Getting the clusters
0159   const auto& clusterHandle = iEvent.getHandle(clustersToken_);
0160   if (!clusterHandle.isValid()) {
0161     edm::LogWarning("Phase2OTValidateCluster") << "No Phase2TrackerCluster1D Collection found in the event. Skipping!";
0162     return;
0163   }
0164   // Getting PixelDigiSimLinks
0165   const auto& pixelSimLinksHandle = iEvent.getHandle(simOTLinksToken_);
0166   if (!pixelSimLinksHandle.isValid()) {
0167     edm::LogWarning("Phase2OTValidateCluster") << "No PixelDigiSimLinks Collection found in the event. Skipping!";
0168     return;
0169   }
0170 
0171   // Number of clusters
0172   std::map<std::string, unsigned int> nPrimarySimHits[3];
0173   std::map<std::string, unsigned int> nOtherSimHits[3];
0174   for (const auto& DSVItr : *clusterHandle) {
0175     // Getting the id of detector unit
0176     uint32_t rawid = DSVItr.detId();
0177     DetId detId(rawid);
0178     const GeomDetUnit* geomDetUnit(tkGeom_->idToDetUnit(detId));
0179     if (!geomDetUnit)
0180       continue;
0181     TrackerGeometry::ModuleType mType = tkGeom_->getDetectorType(detId);
0182 
0183     std::string folderkey = phase2tkutil::getOTHistoId(detId, tTopo_);
0184     for (const auto& clusterItr : DSVItr) {
0185       MeasurementPoint mpCluster(clusterItr.center(), clusterItr.column() + 0.5);
0186       Local3DPoint localPosCluster = geomDetUnit->topology().localPosition(mpCluster);
0187 
0188       // Get simTracks from the cluster
0189       std::vector<unsigned int> clusterSimTrackIds;
0190       for (unsigned int i(0); i < clusterItr.size(); ++i) {
0191         unsigned int channel(Phase2TrackerDigi::pixelToChannel(clusterItr.firstRow() + i, clusterItr.column()));
0192         std::vector<unsigned int> simTrackIds(getSimTrackId(pixelSimLinksHandle, detId, channel));
0193         for (auto it : simTrackIds) {
0194           bool add = true;
0195           for (unsigned int j = 0; j < clusterSimTrackIds.size(); ++j) {
0196             // only save simtrackids that are not present yet
0197             if (it == clusterSimTrackIds.at(j))
0198               add = false;
0199           }
0200           if (add)
0201             clusterSimTrackIds.push_back(it);
0202         }
0203       }
0204       std::sort(clusterSimTrackIds.begin(), clusterSimTrackIds.end());
0205       const PSimHit* closestSimHit = nullptr;
0206       float mind = 1e4;
0207       // Get the SimHit
0208       for (const auto& psimhitCont : simHits) {
0209         for (const auto& simhitIt : *psimhitCont) {
0210           if (rawid == simhitIt.detUnitId()) {
0211             auto it = std::lower_bound(clusterSimTrackIds.begin(), clusterSimTrackIds.end(), simhitIt.trackId());
0212             if (it != clusterSimTrackIds.end() && *it == simhitIt.trackId()) {
0213               float dx = simhitIt.localPosition().x() - localPosCluster.x();
0214               float dy = simhitIt.localPosition().y() - localPosCluster.y();
0215               float dist = dx * dx + dy * dy;
0216               if (!closestSimHit || dist < mind) {
0217                 mind = dist;
0218                 closestSimHit = &simhitIt;
0219               }
0220             }
0221           }
0222         }  //end loop over PSimhitcontainers
0223       }  //end loop over simHits
0224 
0225       if (!closestSimHit)
0226         continue;
0227       // only look at simhits from highpT tracks
0228       auto simTrackIt(simTracks.find(closestSimHit->trackId()));
0229       if (simTrackIt == simTracks.end())
0230         continue;
0231 
0232       Local3DPoint localPosSimHit(closestSimHit->localPosition());
0233       const float deltaX = localPosCluster.x() - localPosSimHit.x();
0234       const float deltaY = localPosCluster.y() - localPosSimHit.y();
0235 
0236       auto layerMEit = layerMEs_.find(folderkey);
0237       if (layerMEit == layerMEs_.end())
0238         continue;
0239 
0240       ClusterMEs& local_mes = layerMEit->second;
0241       if (mType == TrackerGeometry::ModuleType::Ph2PSP) {
0242         local_mes.deltaX_P->Fill(phase2tkutil::cmtomicron * deltaX);
0243         local_mes.deltaY_P->Fill(phase2tkutil::cmtomicron * deltaY);
0244       } else if (mType == TrackerGeometry::ModuleType::Ph2PSS || mType == TrackerGeometry::ModuleType::Ph2SS) {
0245         local_mes.deltaX_S->Fill(phase2tkutil::cmtomicron * deltaX);
0246         local_mes.deltaY_S->Fill(deltaY);
0247       }
0248       // Primary particles only
0249       if (phase2tkutil::isPrimary(simTrackIt->second, closestSimHit)) {
0250         if (mType == TrackerGeometry::ModuleType::Ph2PSP) {
0251           local_mes.deltaX_P_primary->Fill(phase2tkutil::cmtomicron * deltaX);
0252           local_mes.deltaY_P_primary->Fill(phase2tkutil::cmtomicron * deltaY);
0253         } else if (mType == TrackerGeometry::ModuleType::Ph2PSS || mType == TrackerGeometry::ModuleType::Ph2SS) {
0254           local_mes.deltaX_S_primary->Fill(phase2tkutil::cmtomicron * deltaX);
0255           local_mes.deltaY_S_primary->Fill(deltaY);
0256         }
0257       }
0258     }
0259   }
0260 }
0261 
0262 //
0263 // -- Book Histograms
0264 //
0265 void Phase2OTValidateCluster::bookHistograms(DQMStore::IBooker& ibooker,
0266                                              edm::Run const& iRun,
0267                                              edm::EventSetup const& iSetup) {
0268   std::string top_folder = config_.getParameter<std::string>("TopFolderName");
0269   edm::LogInfo("Phase2OTValidateCluster") << " Booking Histograms in: " << top_folder;
0270 
0271   edm::ESWatcher<TrackerDigiGeometryRecord> theTkDigiGeomWatcher;
0272   if (theTkDigiGeomWatcher.check(iSetup)) {
0273     for (auto const& det_u : tkGeom_->detUnits()) {
0274       //Always check TrackerNumberingBuilder before changing this part
0275       if ((det_u->subDetector() == GeomDetEnumerators::SubDetector::P2PXB ||
0276            det_u->subDetector() == GeomDetEnumerators::SubDetector::P2PXEC))
0277         continue;  //continue if Pixel
0278       uint32_t detId_raw = det_u->geographicalId().rawId();
0279       bookLayerHistos(ibooker, detId_raw, top_folder);
0280     }
0281   }
0282 }
0283 
0284 //////////////////Layer Histo/////////////////////////////////
0285 void Phase2OTValidateCluster::bookLayerHistos(DQMStore::IBooker& ibooker, uint32_t det_id, const std::string& subdir) {
0286   std::string folderName = phase2tkutil::getOTHistoId(det_id, tTopo_);
0287   if (folderName.empty()) {
0288     edm::LogWarning("Phase2OTValidateCluster") << ">>>> Invalid histo_id ";
0289     return;
0290   }
0291 
0292   if (layerMEs_.find(folderName) == layerMEs_.end()) {
0293     ibooker.cd();
0294     edm::LogInfo("Phase2TrackerValidateDigi") << " Booking Histograms in: " << subdir + '/' + folderName;
0295     ClusterMEs local_mes;
0296     if (tkGeom_->getDetectorType(det_id) == TrackerGeometry::ModuleType::Ph2PSP) {
0297       ibooker.setCurrentFolder(subdir + '/' + folderName);
0298 
0299       local_mes.deltaX_P =
0300           phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_X_Pixel"), ibooker);
0301 
0302       local_mes.deltaY_P =
0303           phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_Y_Pixel"), ibooker);
0304 
0305       // Puting primary digis in a subfolder
0306       ibooker.setCurrentFolder(subdir + '/' + folderName + "/PrimarySimHits");
0307 
0308       local_mes.deltaX_P_primary =
0309           phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_X_Pixel_Primary"), ibooker);
0310 
0311       local_mes.deltaY_P_primary =
0312           phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_Y_Pixel_Primary"), ibooker);
0313     }
0314     ibooker.setCurrentFolder(subdir + '/' + folderName);
0315 
0316     local_mes.deltaX_S =
0317         phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_X_Strip"), ibooker);
0318 
0319     local_mes.deltaY_S =
0320         phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_Y_Strip"), ibooker);
0321 
0322     // Puting primary digis in a subfolder
0323     ibooker.setCurrentFolder(subdir + '/' + folderName + "/PrimarySimHits");
0324 
0325     local_mes.deltaX_S_primary =
0326         phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_X_Strip_Primary"), ibooker);
0327 
0328     local_mes.deltaY_S_primary =
0329         phase2tkutil::book1DFromPSet(config_.getParameter<edm::ParameterSet>("Delta_Y_Strip_Primary"), ibooker);
0330 
0331     layerMEs_.emplace(folderName, local_mes);
0332   }
0333 }
0334 
0335 std::vector<unsigned int> Phase2OTValidateCluster::getSimTrackId(
0336     const edm::Handle<edm::DetSetVector<PixelDigiSimLink>>& pixelSimLinks, const DetId& detId, unsigned int channel) {
0337   std::vector<unsigned int> retvec;
0338   edm::DetSetVector<PixelDigiSimLink>::const_iterator DSViter(pixelSimLinks->find(detId));
0339   if (DSViter == pixelSimLinks->end())
0340     return retvec;
0341   for (edm::DetSet<PixelDigiSimLink>::const_iterator it = DSViter->data.begin(); it != DSViter->data.end(); ++it) {
0342     if (channel == it->channel()) {
0343       retvec.push_back(it->SimTrackId());
0344     }
0345   }
0346   return retvec;
0347 }
0348 
0349 void Phase2OTValidateCluster::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0350   edm::ParameterSetDescription desc;
0351   //for macro-pixel sensors
0352   std::string mptag = "macro-pixel sensor";
0353   std::string striptag = "strip sensor";
0354   {
0355     edm::ParameterSetDescription psd0;
0356     psd0.add<std::string>("name", "Delta_X_Pixel");
0357     psd0.add<std::string>("title", "#Delta X " + mptag + ";Cluster resolution X coordinate [#mum]");
0358     psd0.add<bool>("switch", true);
0359     psd0.add<double>("xmax", 250);
0360     psd0.add<double>("xmin", -250);
0361     psd0.add<int>("NxBins", 100);
0362     desc.add<edm::ParameterSetDescription>("Delta_X_Pixel", psd0);
0363   }
0364   {
0365     edm::ParameterSetDescription psd0;
0366     psd0.add<std::string>("name", "Delta_Y_Pixel");
0367     psd0.add<std::string>("title", "#Delta Y " + mptag + ";Cluster resolution Y coordinate [#mum]");
0368     psd0.add<bool>("switch", true);
0369     psd0.add<double>("xmin", -1500);
0370     psd0.add<double>("xmax", 1500);
0371     psd0.add<int>("NxBins", 100);
0372     desc.add<edm::ParameterSetDescription>("Delta_Y_Pixel", psd0);
0373   }
0374   {
0375     edm::ParameterSetDescription psd0;
0376     psd0.add<std::string>("name", "Delta_X_Pixel_Primary");
0377     psd0.add<std::string>("title", "#Delta X " + mptag + ";cluster resolution X coordinate [#mum]");
0378     psd0.add<bool>("switch", true);
0379     psd0.add<double>("xmin", -250);
0380     psd0.add<double>("xmax", 250);
0381     psd0.add<int>("NxBins", 100);
0382     desc.add<edm::ParameterSetDescription>("Delta_X_Pixel_Primary", psd0);
0383   }
0384   {
0385     edm::ParameterSetDescription psd0;
0386     psd0.add<std::string>("name", "Delta_Y_Pixel_Primary");
0387     psd0.add<std::string>("title", "#Delta Y " + mptag + ";cluster resolution Y coordinate [#mum]");
0388     psd0.add<bool>("switch", true);
0389     psd0.add<double>("xmin", -1500);
0390     psd0.add<double>("xmax", 1500);
0391     psd0.add<int>("NxBins", 100);
0392     desc.add<edm::ParameterSetDescription>("Delta_Y_Pixel_Primary", psd0);
0393   }
0394 
0395   //strip sensors
0396   {
0397     edm::ParameterSetDescription psd0;
0398     psd0.add<std::string>("name", "Delta_X_Strip");
0399     psd0.add<std::string>("title", "#Delta X " + striptag + ";Cluster resolution X coordinate [#mum]");
0400     psd0.add<bool>("switch", true);
0401     psd0.add<double>("xmin", -250);
0402     psd0.add<double>("xmax", 250);
0403     psd0.add<int>("NxBins", 100);
0404     desc.add<edm::ParameterSetDescription>("Delta_X_Strip", psd0);
0405   }
0406   {
0407     edm::ParameterSetDescription psd0;
0408     psd0.add<std::string>("name", "Delta_Y_Strip");
0409     psd0.add<std::string>("title", "#Delta Y " + striptag + ";Cluster resolution Y coordinate [cm]");
0410     psd0.add<double>("xmin", -5.0);
0411     psd0.add<bool>("switch", true);
0412     psd0.add<double>("xmax", 5.0);
0413     psd0.add<int>("NxBins", 100);
0414     desc.add<edm::ParameterSetDescription>("Delta_Y_Strip", psd0);
0415   }
0416   {
0417     edm::ParameterSetDescription psd0;
0418     psd0.add<std::string>("name", "Delta_X_Strip_Primary");
0419     psd0.add<std::string>("title", "#Delta X " + striptag + ";Cluster resolution X coordinate [#mum]");
0420     psd0.add<bool>("switch", true);
0421     psd0.add<double>("xmin", -250);
0422     psd0.add<double>("xmax", 250);
0423     psd0.add<int>("NxBins", 100);
0424     desc.add<edm::ParameterSetDescription>("Delta_X_Strip_Primary", psd0);
0425   }
0426   {
0427     edm::ParameterSetDescription psd0;
0428     psd0.add<std::string>("name", "Delta_Y_Strip_Primary");
0429     psd0.add<std::string>("title", "#Delta Y " + striptag + ";Cluster resolution Y coordinate [cm]");
0430     psd0.add<double>("xmin", -5.0);
0431     psd0.add<bool>("switch", true);
0432     psd0.add<double>("xmax", 5.0);
0433     psd0.add<int>("NxBins", 100);
0434     desc.add<edm::ParameterSetDescription>("Delta_Y_Strip_Primary", psd0);
0435   }
0436   desc.add<std::string>("TopFolderName", "TrackerPhase2OTClusterV");
0437   desc.add<edm::InputTag>("ClusterSource", edm::InputTag("siPhase2Clusters"));
0438   desc.add<edm::InputTag>("OuterTrackerDigiSimLinkSource", edm::InputTag("simSiPixelDigis", "Tracker"));
0439   desc.add<edm::InputTag>("simtracks", edm::InputTag("g4SimHits"));
0440   desc.add<double>("SimTrackMinPt", 0.0);
0441   desc.add<std::vector<edm::InputTag>>("PSimHitSource",
0442                                        {
0443                                            edm::InputTag("g4SimHits:TrackerHitsTIBLowTof"),
0444                                            edm::InputTag("g4SimHits:TrackerHitsTIBHighTof"),
0445                                            edm::InputTag("g4SimHits:TrackerHitsTIDLowTof"),
0446                                            edm::InputTag("g4SimHits:TrackerHitsTIDHighTof"),
0447                                            edm::InputTag("g4SimHits:TrackerHitsTOBLowTof"),
0448                                            edm::InputTag("g4SimHits:TrackerHitsTOBHighTof"),
0449                                            edm::InputTag("g4SimHits:TrackerHitsTECLowTof"),
0450                                            edm::InputTag("g4SimHits:TrackerHitsTECHighTof"),
0451                                            edm::InputTag("g4SimHits:TrackerHitsPixelBarrelLowTof"),
0452                                            edm::InputTag("g4SimHits:TrackerHitsPixelBarrelHighTof"),
0453                                            edm::InputTag("g4SimHits:TrackerHitsPixelEndcapLowTof"),
0454                                            edm::InputTag("g4SimHits:TrackerHitsPixelEndcapHighTof"),
0455                                        });
0456   descriptions.add("Phase2OTValidateCluster", desc);
0457 }
0458 DEFINE_FWK_MODULE(Phase2OTValidateCluster);