File indexing completed on 2023-03-17 11:19:52
0001 #include "RecoLocalTracker/SubCollectionProducers/interface/ClusterSummaryProducer.h"
0002
0003 ClusterSummaryProducer::ClusterSummaryProducer(const edm::ParameterSet& iConfig)
0004 : doStrips(iConfig.getParameter<bool>("doStrips")),
0005 doPixels(iConfig.getParameter<bool>("doPixels")),
0006 verbose(iConfig.getParameter<bool>("verbose")) {
0007 pixelClusters_ =
0008 consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<edm::InputTag>("pixelClusters"));
0009 stripClusters_ =
0010 consumes<edmNew::DetSetVector<SiStripCluster> >(iConfig.getParameter<edm::InputTag>("stripClusters"));
0011
0012 ClusterSummary::CMSTracker maxEnum = ClusterSummary::STRIP;
0013
0014 std::vector<std::string> wantedsubdets = iConfig.getParameter<std::vector<std::string> >("wantedSubDets");
0015 for (const auto& iS : wantedsubdets) {
0016 ClusterSummary::CMSTracker subdet = ClusterSummary::NVALIDENUMS;
0017 for (int iN = 0; iN < ClusterSummary::NVALIDENUMS; ++iN)
0018 if (ClusterSummary::subDetNames[iN] == iS)
0019 subdet = ClusterSummary::CMSTracker(iN);
0020 if (subdet == ClusterSummary::NVALIDENUMS)
0021 throw cms::Exception("No standard selection: ") << iS;
0022
0023 selectors.push_back(ModuleSelection(DetIdSelector(ClusterSummary::subDetSelections[subdet]), subdet));
0024 if (subdet > maxEnum)
0025 maxEnum = subdet;
0026 if (verbose)
0027 moduleNames.push_back(ClusterSummary::subDetNames[subdet]);
0028 }
0029
0030 std::vector<edm::ParameterSet> wantedusersubdets_ps =
0031 iConfig.getParameter<std::vector<edm::ParameterSet> >("wantedUserSubDets");
0032 for (const auto& iS : wantedusersubdets_ps) {
0033 ClusterSummary::CMSTracker subdet = (ClusterSummary::CMSTracker)iS.getParameter<unsigned int>("detSelection");
0034 std::string detname = iS.getParameter<std::string>("detLabel");
0035 std::vector<std::string> selection = iS.getParameter<std::vector<std::string> >("selection");
0036
0037 if (subdet <= ClusterSummary::NVALIDENUMS)
0038 throw cms::Exception("Already predefined selection: ") << subdet;
0039 if (subdet >= ClusterSummary::NTRACKERENUMS)
0040 throw cms::Exception("Selection is out of range: ") << subdet;
0041
0042 selectors.push_back(ModuleSelection(DetIdSelector(selection), subdet));
0043 if (subdet > maxEnum)
0044 maxEnum = subdet;
0045 if (verbose)
0046 moduleNames.push_back(detname);
0047 }
0048
0049 cCluster = ClusterSummary(maxEnum + 1);
0050 produces<ClusterSummary>().setBranchAlias("trackerClusterSummary");
0051 }
0052
0053 void ClusterSummaryProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0054 using namespace edm;
0055 cCluster.reset();
0056 std::vector<bool> selectedVector(selectors.size(), false);
0057
0058 auto getSelections = [&](const uint32_t detid) {
0059 for (unsigned int iS = 0; iS < selectors.size(); ++iS)
0060 selectedVector[iS] = selectors[iS].first.isSelected(detid);
0061 };
0062 auto fillSelections = [&](const int clusterSize, const float clusterCharge) {
0063 for (unsigned int iS = 0; iS < selectors.size(); ++iS) {
0064 if (!selectedVector[iS])
0065 continue;
0066 const ClusterSummary::CMSTracker module = selectors[iS].second;
0067 cCluster.addNClusByIndex(module, 1);
0068 cCluster.addClusSizeByIndex(module, clusterSize);
0069 cCluster.addClusChargeByIndex(module, clusterCharge);
0070 }
0071 };
0072
0073
0074
0075
0076 if (doStrips) {
0077 edm::Handle<edmNew::DetSetVector<SiStripCluster> > stripClusters;
0078 iEvent.getByToken(stripClusters_, stripClusters);
0079 edmNew::DetSetVector<SiStripCluster>::const_iterator itClusters = stripClusters->begin();
0080 for (; itClusters != stripClusters->end(); ++itClusters) {
0081 getSelections(itClusters->id());
0082 for (edmNew::DetSet<SiStripCluster>::const_iterator cluster = itClusters->begin(); cluster != itClusters->end();
0083 ++cluster) {
0084 const ClusterVariables Summaryinfo(*cluster);
0085 fillSelections(Summaryinfo.clusterSize(), Summaryinfo.charge());
0086 }
0087 }
0088 }
0089
0090
0091
0092
0093 if (doPixels) {
0094 edm::Handle<edmNew::DetSetVector<SiPixelCluster> > pixelClusters;
0095 iEvent.getByToken(pixelClusters_, pixelClusters);
0096 edmNew::DetSetVector<SiPixelCluster>::const_iterator itClusters = pixelClusters->begin();
0097 for (; itClusters != pixelClusters->end(); ++itClusters) {
0098 getSelections(itClusters->id());
0099 for (edmNew::DetSet<SiPixelCluster>::const_iterator cluster = itClusters->begin(); cluster != itClusters->end();
0100 ++cluster) {
0101 fillSelections(cluster->size(), float(cluster->charge()) / 1000.);
0102 }
0103 }
0104 }
0105
0106
0107
0108
0109 if (verbose) {
0110 for (const auto& iS : selectors) {
0111 const ClusterSummary::CMSTracker module = iS.second;
0112 if (cCluster.getNClusByIndex(module) != 0) {
0113 edm::LogInfo("ClusterSummaryProducer")
0114 << "n" << moduleNames[module] << ", avg size, avg charge = " << cCluster.getNClusByIndex(module) << ", "
0115 << cCluster.getClusSizeByIndex(module) / cCluster.getNClusByIndex(module) << ", "
0116 << cCluster.getClusChargeByIndex(module) / cCluster.getNClusByIndex(module) << std::endl;
0117 } else {
0118 edm::LogInfo("ClusterSummaryProducer")
0119 << "n" << moduleNames[module] << ", avg size, avg charge = " << cCluster.getNClusByIndex(module) << ", "
0120 << "0"
0121 << ", "
0122 << "0" << std::endl;
0123 }
0124 edm::LogInfo("ClusterSummaryProducer") << "-------------------------------------------------------" << std::endl;
0125 }
0126 }
0127
0128
0129 auto result = std::make_unique<ClusterSummary>();
0130
0131 result->copyNonEmpty(cCluster);
0132 iEvent.put(std::move(result));
0133 }
0134
0135 void ClusterSummaryProducer::beginStream(edm::StreamID) {
0136 if (!verbose)
0137 return;
0138 edm::LogInfo("ClusterSummaryProducer") << "+++++++++++++++++++++++++++++++ " << std::endl << "Getting info on ";
0139 for (const auto& iS : moduleNames) {
0140 edm::LogInfo("ClusterSummaryProducer") << iS << " ";
0141 }
0142 edm::LogInfo("ClusterSummaryProducer") << std::endl;
0143 }
0144
0145 #include "FWCore/PluginManager/interface/ModuleDef.h"
0146 #include "FWCore/Framework/interface/MakerMacros.h"
0147 DEFINE_FWK_MODULE(ClusterSummaryProducer);