SiStripMonitorApproximateCluster/plugins/SiStripMonitorApproximateCluster.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    DQM/SiStripMonitorApproximateCluster
0004 // Class:      SiStripMonitorApproximateCluster
0005 //
0006 /**\class SiStripMonitorApproximateCluster SiStripMonitorApproximateCluster.cc DQM/SiStripMonitorApproximateCluster/plugins/SiStripMonitorApproximateCluster.cc
0007
0008  Description: Monitor SiStripApproximateClusters and on-demand compare properties with original SiStripClusters
0009
0010 */
0011 //
0012 // Original Author:  Marco Musich
0013 //         Created:  Thu, 08 Dec 2022 20:51:10 GMT
0014 //
0015 //
0016
0017 #include <string>
0018 // for string manipulations
0019 #include <fmt/printf.h>
0020
0021 // user include files
0022 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0023 #include "DQMServices/Core/interface/MonitorElement.h"
0024 #include "DataFormats/Common/interface/DetSet.h"
0025 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0026 #include "DataFormats/DetId/interface/DetIdVector.h"
0027 #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h"
0028 #include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollection.h"
0029 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0030 #include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h"
0031 #include "FWCore/Framework/interface/Event.h"
0032 #include "FWCore/Framework/interface/Frameworkfwd.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0035 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0036 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0037 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0038
0039 //
0040 // class declaration
0041 //
0042
0043 namespace siStripRawPrime {
0044   struct monitorApproxCluster {
0045   public:
0046     monitorApproxCluster()
0047         : h_barycenter_{nullptr}, h_width_{nullptr}, h_avgCharge_{nullptr}, h_isSaturated_{nullptr}, isBooked_{false} {}
0048
0049     void fill(const SiStripApproximateCluster& cluster) {
0050       h_barycenter_->Fill(cluster.barycenter());
0051       h_width_->Fill(cluster.width());
0052       h_charge_->Fill(cluster.width() * cluster.avgCharge());  // see SiStripCluster constructor
0053       h_avgCharge_->Fill(cluster.avgCharge());
0054       h_isSaturated_->Fill(cluster.isSaturated() ? 1 : -1);
0055       h_passFilter_->Fill(cluster.filter() ? 1 : -1);
0056       h_passPeakFilter_->Fill(cluster.peakFilter() ? 1 : -1);
0057     }
0058
0059     void book(dqm::implementation::DQMStore::IBooker& ibook, const std::string& folder) {
0060       ibook.setCurrentFolder(folder);
0061       h_barycenter_ =
0062           ibook.book1D("clusterBarycenter", "cluster barycenter;cluster barycenter;#clusters", 7680., 0., 7680.);
0063       h_width_ = ibook.book1D("clusterWidth", "cluster width;cluster width;#clusters", 128, -0.5, 127.5);
0064       h_avgCharge_ = ibook.book1D(
0065           "clusterAvgCharge", "average strip charge;average strip charge [ADC counts];#clusters", 256, -0.5, 255.5);
0066
0067       h_charge_ = ibook.book1D(
0068           "clusterCharge", "total cluster charge;total cluster charge [ADC counts];#clusters", 300, -0.5, 2999.5);
0069
0070       h_isSaturated_ = ibook.book1D("clusterSaturation", "cluster saturation;is saturated?;#clusters", 3, -1.5, 1.5);
0071       h_isSaturated_->setBinLabel(1, "Not saturated");
0072       h_isSaturated_->setBinLabel(3, "Saturated");
0073
0074       h_passFilter_ = ibook.book1D("filter", "filter;passes filter?;#clusters", 3, -1.5, 1.5);
0075       h_passFilter_->setBinLabel(1, "No");
0076       h_passFilter_->setBinLabel(3, "Yes");
0077
0078       h_passPeakFilter_ = ibook.book1D("peakFilter", "peak filter;passes peak filter?;#clusters", 3, -1.5, 1.5);
0079       h_passPeakFilter_->setBinLabel(1, "No");
0080       h_passPeakFilter_->setBinLabel(3, "Yes");
0081
0082       isBooked_ = true;
0083     }
0084
0085     bool isBooked() { return isBooked_; }
0086
0087   private:
0088     // approximate clusters
0089     dqm::reco::MonitorElement* h_barycenter_;
0090     dqm::reco::MonitorElement* h_width_;
0091     dqm::reco::MonitorElement* h_avgCharge_;
0092     dqm::reco::MonitorElement* h_charge_;
0093     dqm::reco::MonitorElement* h_isSaturated_;
0094     dqm::reco::MonitorElement* h_passFilter_;
0095     dqm::reco::MonitorElement* h_passPeakFilter_;
0096
0097     bool isBooked_;
0098   };
0099 }  // namespace siStripRawPrime
0100
0101 class SiStripMonitorApproximateCluster : public DQMEDAnalyzer {
0102 public:
0103   explicit SiStripMonitorApproximateCluster(const edm::ParameterSet&);
0104   ~SiStripMonitorApproximateCluster() override = default;
0105
0106   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0107
0108 private:
0109   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0110
0111   void analyze(const edm::Event&, const edm::EventSetup&) override;
0112
0113   // ------------ member data ------------
0114   std::string folder_;
0115   bool compareClusters_;
0116   MonitorElement* h_nclusters_;
0117
0118   siStripRawPrime::monitorApproxCluster allClusters{};
0119   siStripRawPrime::monitorApproxCluster matchedClusters{};
0120   siStripRawPrime::monitorApproxCluster unMatchedClusters{};
0121
0122   // FED errors
0123   dqm::reco::MonitorElement* h_numberFEDErrors_;
0124
0125   // for comparisons
0126   MonitorElement* h_isMatched_{nullptr};
0127   MonitorElement* h_deltaBarycenter_{nullptr};
0128   MonitorElement* h_deltaSize_{nullptr};
0129   MonitorElement* h_deltaCharge_{nullptr};
0130   MonitorElement* h_deltaFirstStrip_{nullptr};
0131   MonitorElement* h_deltaEndStrip_{nullptr};
0132
0133   MonitorElement* h2_CompareBarycenter_{nullptr};
0134   MonitorElement* h2_CompareSize_{nullptr};
0135   MonitorElement* h2_CompareCharge_{nullptr};
0136   MonitorElement* h2_CompareFirstStrip_{nullptr};
0137   MonitorElement* h2_CompareEndStrip_{nullptr};
0138
0139   // Event Data
0140   const edm::EDGetTokenT<SiStripApproximateClusterCollection> approxClustersToken_;
0141   const edm::EDGetTokenT<DetIdVector> stripFEDErrorsToken_;
0142
0143   edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster>> stripClustersToken_;
0144   const edmNew::DetSetVector<SiStripCluster>* stripClusterCollection_;
0145
0146   // Event Setup Data
0147   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_;
0148 };
0149
0150 //
0151 // constructors and destructor
0152 //
0153 SiStripMonitorApproximateCluster::SiStripMonitorApproximateCluster(const edm::ParameterSet& iConfig)
0154     : folder_(iConfig.getParameter<std::string>("folder")),
0155       compareClusters_(iConfig.getParameter<bool>("compareClusters")),
0156       // Poducer name of input StripClusterCollection
0157       approxClustersToken_(
0158           consumes<SiStripApproximateClusterCollection>(iConfig.getParameter<edm::InputTag>("ApproxClustersProducer"))),
0159       stripFEDErrorsToken_(consumes<DetIdVector>(iConfig.getParameter<edm::InputTag>("SiStripFEDErrorVector"))) {
0160   tkGeomToken_ = esConsumes();
0161   if (compareClusters_) {
0162     stripClustersToken_ =
0163         consumes<edmNew::DetSetVector<SiStripCluster>>(iConfig.getParameter<edm::InputTag>("ClustersProducer"));
0164   }
0165   stripClusterCollection_ = nullptr;
0166 }
0167
0168 //
0169 // member functions
0170 //
0171
0172 // ------------ method called for each event  ------------
0173 void SiStripMonitorApproximateCluster::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0174   using namespace edm;
0175
0176   const auto& tkGeom = &iSetup.getData(tkGeomToken_);
0177   const auto tkDets = tkGeom->dets();
0178
0179   // get SiStripApproximateClusterCollection from Event
0180   edm::Handle<SiStripApproximateClusterCollection> approx_cluster_detsetvector;
0181   iEvent.getByToken(approxClustersToken_, approx_cluster_detsetvector);
0182   if (!approx_cluster_detsetvector.isValid()) {
0183     edm::LogError("SiStripMonitorApproximateCluster")
0184         << "SiStripApproximate cluster collection is not valid!" << std::endl;
0185
0186     // if approximate clusters collection not available, then early return
0187     return;
0188   }
0189
0190   // get the DetIdVector of the SiStrip FED Errors
0191   edm::Handle<DetIdVector> sistripFEDErrorsVector;
0192   iEvent.getByToken(stripFEDErrorsToken_, sistripFEDErrorsVector);
0193   if (!sistripFEDErrorsVector.isValid()) {
0194     edm::LogError("SiStripMonitorApproximateCluster")
0195         << " DetIdVector collection of SiStrip FED errors is not valid!" << std::endl;
0196
0197     // if approximate clusters collection not available, then early return
0198     return;
0199   }
0200
0201   // if requested to perform the comparison
0202   if (compareClusters_) {
0203     // get collection of DetSetVector of clusters from Event
0204     edm::Handle<edmNew::DetSetVector<SiStripCluster>> cluster_detsetvector;
0205     iEvent.getByToken(stripClustersToken_, cluster_detsetvector);
0206     if (!cluster_detsetvector.isValid()) {
0207       edm::LogError("SiStripMonitorApproximateCluster")
0208           << "Requested to perform comparison, but regular SiStrip cluster collection is not valid!" << std::endl;
0209       return;
0210     } else {
0211       stripClusterCollection_ = cluster_detsetvector.product();
0212     }
0213   }
0214
0215   int nApproxClusters{0};
0216   const SiStripApproximateClusterCollection* clusterCollection = approx_cluster_detsetvector.product();
0217
0218   for (const auto& detClusters : *clusterCollection) {
0219     edmNew::DetSet<SiStripCluster> strip_clusters_detset;
0220     const auto& detid = detClusters.id();  // get the detid of the current detset
0221
0222     // starts here comaparison with regular clusters
0223     if (compareClusters_) {
0224       if (stripClusterCollection_->empty()) {
0225         edm::LogWarning("SiStripMonitorApproximateCluster")
0226             << "Input SiStrip Cluster collecction was empty, skipping event! " << std::endl;
0227         return;
0228       }
0229
0230       edmNew::DetSetVector<SiStripCluster>::const_iterator isearch =
0231           stripClusterCollection_->find(detid);  // search clusters of same detid
0232
0233       // protect against a missing match
0234       if (isearch != stripClusterCollection_->end()) {
0235         strip_clusters_detset = (*isearch);
0236       } else {
0237         edm::LogWarning("SiStripMonitorApproximateCluster")
0238             << "No edmNew::DetSet<SiStripCluster> was found for detid " << detid << std::endl;
0239       }
0240     }
0241
0242     for (const auto& cluster : detClusters) {
0243       nApproxClusters++;
0244
0245       // fill the full cluster collection
0246       allClusters.fill(cluster);
0247
0248       if (compareClusters_ && !strip_clusters_detset.empty()) {
0249         // build the converted cluster for the matching
0250         uint16_t nStrips{0};
0251         auto det = std::find_if(tkDets.begin(), tkDets.end(), [detid](auto& elem) -> bool {
0252           return (elem->geographicalId().rawId() == detid);
0253         });
0254         const StripTopology& p = dynamic_cast<const StripGeomDetUnit*>(*det)->specificTopology();
0255         nStrips = p.nstrips() - 1;
0256
0257         const auto convertedCluster = SiStripCluster(cluster, nStrips);
0258
0259         float distance{9999.};
0260         const SiStripCluster* closestCluster{nullptr};
0261         for (const auto& stripCluster : strip_clusters_detset) {
0262           // by construction the approximated cluster width has same
0263           // size as the original cluster
0264
0265           if (cluster.width() != stripCluster.size()) {
0266             continue;
0267           }
0268
0269           float deltaBarycenter = convertedCluster.barycenter() - stripCluster.barycenter();
0270           if (std::abs(deltaBarycenter) < distance) {
0271             closestCluster = &stripCluster;
0272             distance = std::abs(deltaBarycenter);
0273           }
0274         }
0275
0276         // Matching criteria:
0277         // - if exists a closest cluster in the DetId
0278         // - the size coincides with the original one
0279         if (closestCluster) {
0280           // comparisong plots
0281           h_deltaBarycenter_->Fill(closestCluster->barycenter() - convertedCluster.barycenter());
0282           h_deltaSize_->Fill(closestCluster->size() - convertedCluster.size());
0283           h_deltaCharge_->Fill(closestCluster->charge() - convertedCluster.charge());
0284           h_deltaFirstStrip_->Fill(closestCluster->firstStrip() - convertedCluster.firstStrip());
0285           h_deltaEndStrip_->Fill(closestCluster->endStrip() - convertedCluster.endStrip());
0286
0287           h2_CompareBarycenter_->Fill(closestCluster->barycenter(), convertedCluster.barycenter());
0288           h2_CompareSize_->Fill(closestCluster->size(), convertedCluster.size());
0289           h2_CompareCharge_->Fill(closestCluster->charge(), convertedCluster.charge());
0290           h2_CompareFirstStrip_->Fill(closestCluster->firstStrip(), convertedCluster.firstStrip());
0291           h2_CompareEndStrip_->Fill(closestCluster->endStrip(), convertedCluster.endStrip());
0292
0293           // monitoring plots
0294           matchedClusters.fill(cluster);
0295           h_isMatched_->Fill(1);
0296         } else {
0297           // monitoring plots
0298           unMatchedClusters.fill(cluster);
0299           h_isMatched_->Fill(-1);
0300         }
0301       }  // if we're doing the comparison cluster by cluster
0302
0303     }  // loop on clusters in a detset
0304   }  // loop on the detset vector
0305
0306   h_nclusters_->Fill(nApproxClusters);
0307   h_numberFEDErrors_->Fill((*sistripFEDErrorsVector).size());
0308 }
0309
0310 void SiStripMonitorApproximateCluster::bookHistograms(DQMStore::IBooker& ibook,
0311                                                       edm::Run const& run,
0312                                                       edm::EventSetup const& iSetup) {
0313   ibook.setCurrentFolder(folder_);
0314   h_nclusters_ = ibook.book1D("numberOfClusters", "total N. of clusters;N. of clusters;# events", 500., 0., 500000.);
0315   h_numberFEDErrors_ = ibook.book1D(
0316       "numberOfFEDErrors", "number of SiStrip modules in FED Error;N. of Modules in Error; #events", 100, 0, 1000);
0317
0318   allClusters.book(ibook, folder_);
0319
0320   //  for comparisons
0321   if (compareClusters_) {
0322     // book monitoring for matche and unmatched clusters separately
0323     matchedClusters.book(ibook, fmt::format("{}/MatchedClusters", folder_));
0324     unMatchedClusters.book(ibook, fmt::format("{}/UnmatchedClusters", folder_));
0325
0326     // 1D histograms
0327     ibook.setCurrentFolder(fmt::format("{}/ClusterComparisons", folder_));
0328     h_deltaBarycenter_ =
0329         ibook.book1D("deltaBarycenter", "#Delta barycenter;#Delta barycenter;cluster pairs", 101, -50.5, 50.5);
0330     h_deltaSize_ = ibook.book1D("deltaSize", "#Delta size;#Delta size;cluster pairs", 201, -100.5, 100.5);
0331     h_deltaCharge_ = ibook.book1D("deltaCharge", "#Delta charge;#Delta charge;cluster pairs", 401, -200.5, 200.5);
0332
0333     h_deltaFirstStrip_ =
0334         ibook.book1D("deltaFirstStrip", "#Delta FirstStrip; #Delta firstStrip;cluster pairs", 101, -50.5, 50.5);
0335     h_deltaEndStrip_ =
0336         ibook.book1D("deltaEndStrip", "#Delta EndStrip; #Delta endStrip; cluster pairs", 101, -50.5, 50.5);
0337
0338     // geometric constants
0339     constexpr int maxNStrips = 6 * sistrip::STRIPS_PER_APV;
0340     constexpr float minStrip = -0.5f;
0341     constexpr float maxStrip = maxNStrips - 0.5f;
0342
0343     // cluster constants
0344     constexpr float maxCluSize = 100;
0345     constexpr float maxCluCharge = 3000;
0346
0347     // 2D histograms (use TH2I for counts to limit memory allocation)
0348     std::string toRep = "SiStrip Cluster Barycenter";
0349     h2_CompareBarycenter_ = ibook.book2I("compareBarycenter",
0350                                          fmt::sprintf("; %s;Approx %s", toRep, toRep),
0351                                          maxNStrips,
0352                                          minStrip,
0353                                          maxStrip,
0354                                          maxNStrips,
0355                                          minStrip,
0356                                          maxStrip);
0357
0358     toRep = "SiStrip Cluster Size";
0359     h2_CompareSize_ = ibook.book2I("compareSize",
0360                                    fmt::sprintf("; %s;Approx %s", toRep, toRep),
0361                                    maxCluSize,
0362                                    -0.5f,
0363                                    maxCluSize - 0.5f,
0364                                    maxCluSize,
0365                                    -0.5f,
0366                                    maxCluSize - 0.5f);
0367
0368     toRep = "SiStrip Cluster Charge";
0369     h2_CompareCharge_ = ibook.book2I("compareCharge",
0370                                      fmt::sprintf("; %s;Approx %s", toRep, toRep),
0371                                      (maxCluCharge / 5),
0372                                      -0.5f,
0373                                      maxCluCharge - 0.5f,
0374                                      (maxCluCharge / 5),
0375                                      -0.5f,
0376                                      maxCluCharge - 0.5f);
0377
0378     toRep = "SiStrip Cluster First Strip number";
0379     h2_CompareFirstStrip_ = ibook.book2I("compareFirstStrip",
0380                                          fmt::sprintf("; %s;Approx %s", toRep, toRep),
0381                                          maxNStrips,
0382                                          minStrip,
0383                                          maxStrip,
0384                                          maxNStrips,
0385                                          minStrip,
0386                                          maxStrip);
0387
0388     toRep = "SiStrip Cluster Last Strip number";
0389     h2_CompareEndStrip_ = ibook.book2I("compareLastStrip",
0390                                        fmt::sprintf("; %s;Approx %s", toRep, toRep),
0391                                        maxNStrips,
0392                                        minStrip,
0393                                        maxStrip,
0394                                        maxNStrips,
0395                                        minStrip,
0396                                        maxStrip);
0397
0398     h_isMatched_ = ibook.book1D("isClusterMatched", "cluster matching;is matched?;#clusters", 3, -1.5, 1.5);
0399     h_isMatched_->setBinLabel(1, "Not matched");
0400     h_isMatched_->setBinLabel(3, "Matched");
0401   }
0402 }
0403
0404 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0405 void SiStripMonitorApproximateCluster::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0406   edm::ParameterSetDescription desc;
0407   desc.setComment("Monitor SiStripApproximateCluster collection and compare with regular SiStrip clusters");
0408   desc.add<bool>("compareClusters", false)->setComment("if true, will compare with regualr Strip clusters");
0409   desc.add<edm::InputTag>("ApproxClustersProducer", edm::InputTag("hltSiStripClusters2ApproxClusters"))
0410       ->setComment("approxmate clusters collection");
0411   desc.add<edm::InputTag>("SiStripFEDErrorVector", edm::InputTag("hltSiStripRawToDigi"))
0412       ->setComment("SiStrip FED Errors collection");
0413   desc.add<edm::InputTag>("ClustersProducer", edm::InputTag("hltSiStripClusterizerForRawPrime"))
0414       ->setComment("regular clusters collection");
0415   desc.add<std::string>("folder", "SiStripApproximateClusters")->setComment("Top Level Folder");
0416   descriptions.addWithDefaultLabel(desc);
0417 }
0418
0419 // define this as a plug-in
0420 DEFINE_FWK_MODULE(SiStripMonitorApproximateCluster);