CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiStripMonitorApproximateCluster.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: DQM/SiStripMonitorApproximateCluster
4 // Class: SiStripMonitorApproximateCluster
5 //
11 //
12 // Original Author: Marco Musich
13 // Created: Thu, 08 Dec 2022 20:51:10 GMT
14 //
15 //
16 
17 #include <string>
18 // for string manipulations
19 #include <fmt/printf.h>
20 
21 // user include files
38 
39 //
40 // class declaration
41 //
42 
43 namespace siStripRawPrime {
45  public:
47  : h_barycenter_{nullptr}, h_width_{nullptr}, h_avgCharge_{nullptr}, h_isSaturated_{nullptr}, isBooked_{false} {}
48 
49  void fill(const SiStripApproximateCluster& cluster) {
50  h_barycenter_->Fill(cluster.barycenter());
51  h_width_->Fill(cluster.width());
52  h_charge_->Fill(cluster.width() * cluster.avgCharge()); // see SiStripCluster constructor
53  h_avgCharge_->Fill(cluster.avgCharge());
54  h_isSaturated_->Fill(cluster.isSaturated() ? 1 : -1);
55  h_passFilter_->Fill(cluster.filter() ? 1 : -1);
56  h_passPeakFilter_->Fill(cluster.peakFilter() ? 1 : -1);
57  }
58 
60  ibook.setCurrentFolder(folder);
62  ibook.book1D("clusterBarycenter", "cluster barycenter;cluster barycenter;#clusters", 7680., 0., 7680.);
63  h_width_ = ibook.book1D("clusterWidth", "cluster width;cluster width;#clusters", 128, -0.5, 127.5);
64  h_avgCharge_ = ibook.book1D(
65  "clusterAvgCharge", "average strip charge;average strip charge [ADC counts];#clusters", 256, -0.5, 255.5);
66 
67  h_charge_ = ibook.book1D(
68  "clusterCharge", "total cluster charge;total cluster charge [ADC counts];#clusters", 300, -0.5, 2999.5);
69 
70  h_isSaturated_ = ibook.book1D("clusterSaturation", "cluster saturation;is saturated?;#clusters", 3, -1.5, 1.5);
71  h_isSaturated_->setBinLabel(1, "Not saturated");
72  h_isSaturated_->setBinLabel(3, "Saturated");
73 
74  h_passFilter_ = ibook.book1D("filter", "filter;passes filter?;#clusters", 3, -1.5, 1.5);
75  h_passFilter_->setBinLabel(1, "No");
76  h_passFilter_->setBinLabel(3, "Yes");
77 
78  h_passPeakFilter_ = ibook.book1D("peakFilter", "peak filter;passes peak filter?;#clusters", 3, -1.5, 1.5);
80  h_passPeakFilter_->setBinLabel(3, "Yes");
81 
82  isBooked_ = true;
83  }
84 
85  bool isBooked() { return isBooked_; }
86 
87  private:
88  // approximate clusters
96 
97  bool isBooked_;
98  };
99 } // namespace siStripRawPrime
100 
102 public:
104  ~SiStripMonitorApproximateCluster() override = default;
105 
106  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
107 
108 private:
109  void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
110 
111  void analyze(const edm::Event&, const edm::EventSetup&) override;
112 
113  // ------------ member data ------------
117 
121 
122  // FED errors
124 
125  // for comparisons
132 
138 
139  // Event Data
142 
145 
146  // Event Setup Data
148 };
149 
150 //
151 // constructors and destructor
152 //
154  : folder_(iConfig.getParameter<std::string>("folder")),
155  compareClusters_(iConfig.getParameter<bool>("compareClusters")),
156  // Poducer name of input StripClusterCollection
157  approxClustersToken_(
158  consumes<SiStripApproximateClusterCollection>(iConfig.getParameter<edm::InputTag>("ApproxClustersProducer"))),
159  stripFEDErrorsToken_(consumes<DetIdVector>(iConfig.getParameter<edm::InputTag>("SiStripFEDErrorVector"))) {
161  if (compareClusters_) {
163  consumes<edmNew::DetSetVector<SiStripCluster>>(iConfig.getParameter<edm::InputTag>("ClustersProducer"));
164  }
165  stripClusterCollection_ = nullptr;
166 }
167 
168 //
169 // member functions
170 //
171 
172 // ------------ method called for each event ------------
174  using namespace edm;
175 
176  const auto& tkGeom = &iSetup.getData(tkGeomToken_);
177  const auto tkDets = tkGeom->dets();
178 
179  // get SiStripApproximateClusterCollection from Event
180  edm::Handle<SiStripApproximateClusterCollection> approx_cluster_detsetvector;
181  iEvent.getByToken(approxClustersToken_, approx_cluster_detsetvector);
182  if (!approx_cluster_detsetvector.isValid()) {
183  edm::LogError("SiStripMonitorApproximateCluster")
184  << "SiStripApproximate cluster collection is not valid!" << std::endl;
185 
186  // if approximate clusters collection not available, then early return
187  return;
188  }
189 
190  // get the DetIdVector of the SiStrip FED Errors
191  edm::Handle<DetIdVector> sistripFEDErrorsVector;
192  iEvent.getByToken(stripFEDErrorsToken_, sistripFEDErrorsVector);
193  if (!sistripFEDErrorsVector.isValid()) {
194  edm::LogError("SiStripMonitorApproximateCluster")
195  << " DetIdVector collection of SiStrip FED errors is not valid!" << std::endl;
196 
197  // if approximate clusters collection not available, then early return
198  return;
199  }
200 
201  // if requested to perform the comparison
202  if (compareClusters_) {
203  // get collection of DetSetVector of clusters from Event
205  iEvent.getByToken(stripClustersToken_, cluster_detsetvector);
206  if (!cluster_detsetvector.isValid()) {
207  edm::LogError("SiStripMonitorApproximateCluster")
208  << "Requested to perform comparison, but regular SiStrip cluster collection is not valid!" << std::endl;
209  return;
210  } else {
211  stripClusterCollection_ = cluster_detsetvector.product();
212  }
213  }
214 
215  int nApproxClusters{0};
216  const SiStripApproximateClusterCollection* clusterCollection = approx_cluster_detsetvector.product();
217 
218  for (const auto& detClusters : *clusterCollection) {
219  edmNew::DetSet<SiStripCluster> strip_clusters_detset;
220  const auto& detid = detClusters.id(); // get the detid of the current detset
221 
222  // starts here comaparison with regular clusters
223  if (compareClusters_) {
225  stripClusterCollection_->find(detid); // search clusters of same detid
226 
227  // protect against a missing match
228  if (isearch != stripClusterCollection_->end())
229  strip_clusters_detset = (*isearch);
230  }
231 
232  for (const auto& cluster : detClusters) {
233  nApproxClusters++;
234 
235  // fill the full cluster collection
236  allClusters.fill(cluster);
237 
238  if (compareClusters_ && !strip_clusters_detset.empty()) {
239  // build the converted cluster for the matching
240  uint16_t nStrips{0};
241  auto det = std::find_if(tkDets.begin(), tkDets.end(), [detid](auto& elem) -> bool {
242  return (elem->geographicalId().rawId() == detid);
243  });
244  const StripTopology& p = dynamic_cast<const StripGeomDetUnit*>(*det)->specificTopology();
245  nStrips = p.nstrips() - 1;
246 
247  const auto convertedCluster = SiStripCluster(cluster, nStrips);
248 
249  float distance{9999.};
250  const SiStripCluster* closestCluster{nullptr};
251  for (const auto& stripCluster : strip_clusters_detset) {
252  // by construction the approximated cluster width has same
253  // size as the original cluster
254 
255  if (cluster.width() != stripCluster.size()) {
256  continue;
257  }
258 
259  float deltaBarycenter = convertedCluster.barycenter() - stripCluster.barycenter();
260  if (std::abs(deltaBarycenter) < distance) {
261  closestCluster = &stripCluster;
262  distance = std::abs(deltaBarycenter);
263  }
264  }
265 
266  // Matching criteria:
267  // - if exists a closest cluster in the DetId
268  // - the size coincides with the original one
269  if (closestCluster) {
270  // comparisong plots
271  h_deltaBarycenter_->Fill(closestCluster->barycenter() - convertedCluster.barycenter());
272  h_deltaSize_->Fill(closestCluster->size() - convertedCluster.size());
273  h_deltaCharge_->Fill(closestCluster->charge() - convertedCluster.charge());
274  h_deltaFirstStrip_->Fill(closestCluster->firstStrip() - convertedCluster.firstStrip());
275  h_deltaEndStrip_->Fill(closestCluster->endStrip() - convertedCluster.endStrip());
276 
277  h2_CompareBarycenter_->Fill(closestCluster->barycenter(), convertedCluster.barycenter());
278  h2_CompareSize_->Fill(closestCluster->size(), convertedCluster.size());
279  h2_CompareCharge_->Fill(closestCluster->charge(), convertedCluster.charge());
280  h2_CompareFirstStrip_->Fill(closestCluster->firstStrip(), convertedCluster.firstStrip());
281  h2_CompareEndStrip_->Fill(closestCluster->endStrip(), convertedCluster.endStrip());
282 
283  // monitoring plots
284  matchedClusters.fill(cluster);
285  h_isMatched_->Fill(1);
286  } else {
287  // monitoring plots
288  unMatchedClusters.fill(cluster);
289  h_isMatched_->Fill(-1);
290  }
291  } // if we're doing the comparison cluster by cluster
292 
293  } // loop on clusters in a detset
294  } // loop on the detset vector
295 
296  h_nclusters_->Fill(nApproxClusters);
297  h_numberFEDErrors_->Fill((*sistripFEDErrorsVector).size());
298 }
299 
301  edm::Run const& run,
302  edm::EventSetup const& iSetup) {
303  ibook.setCurrentFolder(folder_);
304  h_nclusters_ = ibook.book1D("numberOfClusters", "total N. of clusters;N. of clusters;# events", 500., 0., 500000.);
305  h_numberFEDErrors_ = ibook.book1D(
306  "numberOfFEDErrors", "number of SiStrip modules in FED Error;N. of Modules in Error; #events", 100, 0, 1000);
307 
308  allClusters.book(ibook, folder_);
309 
310  // for comparisons
311  if (compareClusters_) {
312  // book monitoring for matche and unmatched clusters separately
313  matchedClusters.book(ibook, fmt::format("{}/MatchedClusters", folder_));
314  unMatchedClusters.book(ibook, fmt::format("{}/UnmatchedClusters", folder_));
315 
316  // 1D histograms
317  ibook.setCurrentFolder(fmt::format("{}/ClusterComparisons", folder_));
319  ibook.book1D("deltaBarycenter", "#Delta barycenter;#Delta barycenter;cluster pairs", 101, -50.5, 50.5);
320  h_deltaSize_ = ibook.book1D("deltaSize", "#Delta size;#Delta size;cluster pairs", 201, -100.5, 100.5);
321  h_deltaCharge_ = ibook.book1D("deltaCharge", "#Delta charge;#Delta charge;cluster pairs", 401, -200.5, 200.5);
322 
324  ibook.book1D("deltaFirstStrip", "#Delta FirstStrip; #Delta firstStrip;cluster pairs", 101, -50.5, 50.5);
326  ibook.book1D("deltaEndStrip", "#Delta EndStrip; #Delta endStrip; cluster pairs", 101, -50.5, 50.5);
327 
328  // geometric constants
329  constexpr int maxNStrips = 6 * sistrip::STRIPS_PER_APV;
330  constexpr float minStrip = -0.5f;
331  constexpr float maxStrip = maxNStrips - 0.5f;
332 
333  // cluster constants
334  constexpr float maxCluSize = 100;
335  constexpr float maxCluCharge = 3000;
336 
337  // 2D histograms (use TH2I for counts to limit memory allocation)
338  std::string toRep = "SiStrip Cluster Barycenter";
339  h2_CompareBarycenter_ = ibook.book2I("compareBarycenter",
340  fmt::sprintf("; %s;Approx %s", toRep, toRep),
341  maxNStrips,
342  minStrip,
343  maxStrip,
344  maxNStrips,
345  minStrip,
346  maxStrip);
347 
348  toRep = "SiStrip Cluster Size";
349  h2_CompareSize_ = ibook.book2I("compareSize",
350  fmt::sprintf("; %s;Approx %s", toRep, toRep),
351  maxCluSize,
352  -0.5f,
353  maxCluSize - 0.5f,
354  maxCluSize,
355  -0.5f,
356  maxCluSize - 0.5f);
357 
358  toRep = "SiStrip Cluster Charge";
359  h2_CompareCharge_ = ibook.book2I("compareCharge",
360  fmt::sprintf("; %s;Approx %s", toRep, toRep),
361  (maxCluCharge / 5),
362  -0.5f,
363  maxCluCharge - 0.5f,
364  (maxCluCharge / 5),
365  -0.5f,
366  maxCluCharge - 0.5f);
367 
368  toRep = "SiStrip Cluster First Strip number";
369  h2_CompareFirstStrip_ = ibook.book2I("compareFirstStrip",
370  fmt::sprintf("; %s;Approx %s", toRep, toRep),
371  maxNStrips,
372  minStrip,
373  maxStrip,
374  maxNStrips,
375  minStrip,
376  maxStrip);
377 
378  toRep = "SiStrip Cluster Last Strip number";
379  h2_CompareEndStrip_ = ibook.book2I("compareLastStrip",
380  fmt::sprintf("; %s;Approx %s", toRep, toRep),
381  maxNStrips,
382  minStrip,
383  maxStrip,
384  maxNStrips,
385  minStrip,
386  maxStrip);
387 
388  h_isMatched_ = ibook.book1D("isClusterMatched", "cluster matching;is matched?;#clusters", 3, -1.5, 1.5);
389  h_isMatched_->setBinLabel(1, "Not matched");
390  h_isMatched_->setBinLabel(3, "Matched");
391  }
392 }
393 
394 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
397  desc.setComment("Monitor SiStripApproximateCluster collection and compare with regular SiStrip clusters");
398  desc.add<bool>("compareClusters", false)->setComment("if true, will compare with regualr Strip clusters");
399  desc.add<edm::InputTag>("ApproxClustersProducer", edm::InputTag("hltSiStripClusters2ApproxClusters"))
400  ->setComment("approxmate clusters collection");
401  desc.add<edm::InputTag>("SiStripFEDErrorVector", edm::InputTag("hltSiStripRawToDigi"))
402  ->setComment("SiStrip FED Errors collection");
403  desc.add<edm::InputTag>("ClustersProducer", edm::InputTag("hltSiStripClusterizerForRawPrime"))
404  ->setComment("regular clusters collection");
405  desc.add<std::string>("folder", "SiStripApproximateClusters")->setComment("Top Level Folder");
406  descriptions.addWithDefaultLabel(desc);
407 }
408 
409 // define this as a plug-in
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
id_type id() const
Definition: DetSetNew.h:61
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
siStripRawPrime::monitorApproxCluster unMatchedClusters
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
std::string folder_
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
const edm::EDGetTokenT< DetIdVector > stripFEDErrorsToken_
T const * product() const
Definition: Handle.h:70
Log< level::Error, false > LogError
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
bool empty() const
Definition: DetSetNew.h:67
nStrips
1.2 is to make the matching window safely the two nearest strips 0.35 is the size of an ME0 chamber i...
const_iterator end(bool update=false) const
void Fill(long long x)
void fill(const SiStripApproximateCluster &cluster)
int iEvent
Definition: GenABIO.cc:224
virtual const StripTopology & specificTopology() const
Returns a reference to the strip proxy topology.
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double f[11][100]
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
edm::EDGetTokenT< edmNew::DetSetVector< SiStripCluster > > stripClustersToken_
virtual void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
std::vector< DetId > DetIdVector
Definition: DetIdVector.h:7
SiStripMonitorApproximateCluster(const edm::ParameterSet &)
~SiStripMonitorApproximateCluster() override=default
siStripRawPrime::monitorApproxCluster allClusters
Constants and enumerated types for FED/FEC systems.
const edm::EDGetTokenT< SiStripApproximateClusterCollection > approxClustersToken_
const edmNew::DetSetVector< SiStripCluster > * stripClusterCollection_
bool isValid() const
Definition: HandleBase.h:70
HLT enums.
const_iterator find(id_type i, bool update=false) const
static const uint16_t STRIPS_PER_APV
float barycenter() const
void analyze(const edm::Event &, const edm::EventSetup &) override
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
MonitorElement * book2I(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:296
void book(dqm::implementation::DQMStore::IBooker &ibook, const std::string &folder)
Definition: Run.h:45
siStripRawPrime::monitorApproxCluster matchedClusters