Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-12-24 02:22:28

0001 //// -*- C++ -*-
0002 //
0003 // Package:    HGCalHitValidation
0004 // Class:      HGCalHitValidation
0005 //
0006 /**\class HGCalHitValidation HGCalHitValidation.cc Validation/HGCalValidation/plugins/HGCalHitValidation.cc
0007 
0008  Description: [one line class summary]
0009 
0010  Implementation:
0011     [Notes on implementation]
0012 */
0013 
0014 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0015 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0016 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0017 #include "Geometry/HGCalCommonData/interface/HGCalGeometryMode.h"
0018 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0019 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0020 
0021 #include "DataFormats/Common/interface/Handle.h"
0022 #include "DataFormats/HGCRecHit/interface/HGCRecHit.h"
0023 #include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
0024 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0025 #include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
0026 
0027 #include "FWCore/Framework/interface/Frameworkfwd.h"
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/ESHandle.h"
0030 #include "FWCore/ServiceRegistry/interface/Service.h"
0031 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0032 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0033 #include "FWCore/Utilities/interface/Exception.h"
0034 #include "FWCore/Utilities/interface/EDGetToken.h"
0035 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0036 #include "DQMServices/Core/interface/DQMStore.h"
0037 
0038 #include "SimG4CMS/Calo/interface/HGCNumberingScheme.h"
0039 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0040 #include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
0041 #include "SimDataFormats/CaloTest/interface/HGCalTestNumbering.h"
0042 
0043 #include <cmath>
0044 #include <memory>
0045 #include <iostream>
0046 #include <string>
0047 #include <vector>
0048 
0049 //#define EDM_ML_DEBUG
0050 
0051 class HGCalHitValidation : public DQMEDAnalyzer {
0052 public:
0053   explicit HGCalHitValidation(const edm::ParameterSet&);
0054   ~HGCalHitValidation() override;
0055   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0056 
0057 protected:
0058   typedef std::tuple<float, float, float, float> HGCHitTuple;
0059 
0060   void dqmBeginRun(edm::Run const&, edm::EventSetup const&) override;
0061   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0062   void analyze(const edm::Event&, const edm::EventSetup&) override;
0063   void analyzeHGCalSimHit(edm::Handle<std::vector<PCaloHit>> const& simHits,
0064                           int idet,
0065                           MonitorElement* hist,
0066                           std::map<unsigned int, HGCHitTuple>&);
0067 
0068 private:
0069   //HGC Geometry
0070   std::vector<const HGCalDDDConstants*> hgcCons_;
0071   std::vector<const HGCalGeometry*> hgcGeometry_;
0072   std::vector<std::string> geometrySource_;
0073   std::vector<int> ietaExcludeBH_;
0074 
0075   edm::InputTag eeSimHitSource, fhSimHitSource, bhSimHitSource;
0076   edm::EDGetTokenT<std::vector<PCaloHit>> eeSimHitToken_;
0077   edm::EDGetTokenT<std::vector<PCaloHit>> fhSimHitToken_;
0078   edm::EDGetTokenT<std::vector<PCaloHit>> bhSimHitToken_;
0079   edm::EDGetTokenT<HGCeeRecHitCollection> eeRecHitToken_;
0080   edm::EDGetTokenT<HGChefRecHitCollection> fhRecHitToken_;
0081   edm::EDGetTokenT<HGChebRecHitCollection> bhRecHitToken_;
0082   std::vector<edm::ESGetToken<HGCalDDDConstants, IdealGeometryRecord>> tok_ddd_;
0083   std::vector<edm::ESGetToken<HGCalGeometry, IdealGeometryRecord>> tok_geom_;
0084 
0085   //histogram related stuff
0086   MonitorElement *heedzVsZ, *heedyVsY, *heedxVsX;
0087   MonitorElement *hefdzVsZ, *hefdyVsY, *hefdxVsX;
0088   MonitorElement *hebdzVsZ, *hebdPhiVsPhi, *hebdEtaVsEta;
0089 
0090   MonitorElement *heeRecVsSimZ, *heeRecVsSimY, *heeRecVsSimX;
0091   MonitorElement *hefRecVsSimZ, *hefRecVsSimY, *hefRecVsSimX;
0092   MonitorElement *hebRecVsSimZ, *hebRecVsSimY, *hebRecVsSimX;
0093   MonitorElement *heeEnSimRec, *hefEnSimRec, *hebEnSimRec;
0094 
0095   MonitorElement *hebEnRec, *hebEnSim;
0096   MonitorElement *hefEnRec, *hefEnSim;
0097   MonitorElement *heeEnRec, *heeEnSim;
0098 };
0099 
0100 HGCalHitValidation::HGCalHitValidation(const edm::ParameterSet& cfg) {
0101   geometrySource_ = cfg.getParameter<std::vector<std::string>>("geometrySource");
0102   eeSimHitToken_ = consumes<std::vector<PCaloHit>>(cfg.getParameter<edm::InputTag>("eeSimHitSource"));
0103   fhSimHitToken_ = consumes<std::vector<PCaloHit>>(cfg.getParameter<edm::InputTag>("fhSimHitSource"));
0104   bhSimHitToken_ = consumes<std::vector<PCaloHit>>(cfg.getParameter<edm::InputTag>("bhSimHitSource"));
0105   eeRecHitToken_ = consumes<HGCeeRecHitCollection>(cfg.getParameter<edm::InputTag>("eeRecHitSource"));
0106   fhRecHitToken_ = consumes<HGChefRecHitCollection>(cfg.getParameter<edm::InputTag>("fhRecHitSource"));
0107   bhRecHitToken_ = consumes<HGChebRecHitCollection>(cfg.getParameter<edm::InputTag>("bhRecHitSource"));
0108   ietaExcludeBH_ = cfg.getParameter<std::vector<int>>("ietaExcludeBH");
0109 
0110   for (const auto& name : geometrySource_) {
0111     tok_ddd_.emplace_back(
0112         esConsumes<HGCalDDDConstants, IdealGeometryRecord, edm::Transition::BeginRun>(edm::ESInputTag{"", name}));
0113     tok_geom_.emplace_back(
0114         esConsumes<HGCalGeometry, IdealGeometryRecord, edm::Transition::BeginRun>(edm::ESInputTag{"", name}));
0115   }
0116 
0117 #ifdef EDM_ML_DEBUG
0118   edm::LogInfo("HGCalValid") << "Exclude the following " << ietaExcludeBH_.size() << " ieta values from BH plots";
0119   for (unsigned int k = 0; k < ietaExcludeBH_.size(); ++k)
0120     edm::LogInfo("HGCalValid") << " [" << k << "] " << ietaExcludeBH_[k];
0121 #endif
0122 }
0123 
0124 HGCalHitValidation::~HGCalHitValidation() {}
0125 
0126 void HGCalHitValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0127   edm::ParameterSetDescription desc;
0128   std::vector<std::string> source = {"HGCalEESensitive", "HGCalHESiliconSensitive", "HGCalHEScintillatorSensitive"};
0129   desc.add<std::vector<std::string>>("geometrySource", source);
0130   desc.add<edm::InputTag>("eeSimHitSource", edm::InputTag("g4SimHits", "HGCHitsEE"));
0131   desc.add<edm::InputTag>("fhSimHitSource", edm::InputTag("g4SimHits", "HGCHitsHEfront"));
0132   desc.add<edm::InputTag>("bhSimHitSource", edm::InputTag("g4SimHits", "HGCHitsHEback"));
0133   desc.add<edm::InputTag>("eeRecHitSource", edm::InputTag("HGCalRecHit", "HGCEERecHits"));
0134   desc.add<edm::InputTag>("fhRecHitSource", edm::InputTag("HGCalRecHit", "HGCHEFRecHits"));
0135   desc.add<edm::InputTag>("bhRecHitSource", edm::InputTag("HGCalRecHit", "HGCHEBRecHits"));
0136   std::vector<int> dummy;
0137   desc.add<std::vector<int>>("ietaExcludeBH", dummy);
0138   descriptions.add("hgcalHitValidation", desc);
0139 }
0140 
0141 void HGCalHitValidation::bookHistograms(DQMStore::IBooker& iB, edm::Run const&, edm::EventSetup const&) {
0142   iB.setCurrentFolder("HGCAL/HGCalSimHitsV/HitValidation");
0143 
0144   //initiating histograms
0145   heedzVsZ = iB.book2D("heedzVsZ", "", 720, -360, 360, 100, -0.1, 0.1);
0146   heedyVsY = iB.book2D("heedyVsY", "", 400, -200, 200, 100, -0.02, 0.02);
0147   heedxVsX = iB.book2D("heedxVsX", "", 400, -200, 200, 100, -0.02, 0.02);
0148   heeRecVsSimZ = iB.book2D("heeRecVsSimZ", "", 720, -360, 360, 720, -360, 360);
0149   heeRecVsSimY = iB.book2D("heeRecVsSimY", "", 400, -200, 200, 400, -200, 200);
0150   heeRecVsSimX = iB.book2D("heeRecVsSimX", "", 400, -200, 200, 400, -200, 200);
0151 
0152   hefdzVsZ = iB.book2D("hefdzVsZ", "", 820, -410, 410, 100, -0.1, 0.1);
0153   hefdyVsY = iB.book2D("hefdyVsY", "", 400, -200, 200, 100, -0.02, 0.02);
0154   hefdxVsX = iB.book2D("hefdxVsX", "", 400, -200, 200, 100, -0.02, 0.02);
0155   hefRecVsSimZ = iB.book2D("hefRecVsSimZ", "", 820, -410, 410, 820, -410, 410);
0156   hefRecVsSimY = iB.book2D("hefRecVsSimY", "", 400, -200, 200, 400, -200, 200);
0157   hefRecVsSimX = iB.book2D("hefRecVsSimX", "", 400, -200, 200, 400, -200, 200);
0158 
0159   hebdzVsZ = iB.book2D("hebdzVsZ", "", 1080, -540, 540, 100, -1.0, 1.0);
0160   hebdPhiVsPhi = iB.book2D("hebdPhiVsPhi", "", M_PI * 100, -0.5, M_PI + 0.5, 200, -0.2, 0.2);
0161   hebdEtaVsEta = iB.book2D("hebdEtaVsEta", "", 1000, -5, 5, 200, -0.1, 0.1);
0162   hebRecVsSimZ = iB.book2D("hebRecVsSimZ", "", 1080, -540, 540, 1080, -540, 540);
0163   hebRecVsSimY = iB.book2D("hebRecVsSimY", "", 400, -200, 200, 400, -200, 200);
0164   hebRecVsSimX = iB.book2D("hebRecVsSimX", "", 400, -200, 200, 400, -200, 200);
0165 
0166   heeEnRec = iB.book1D("heeEnRec", "", 1000, 0, 10);
0167   heeEnSim = iB.book1D("heeEnSim", "", 1000, 0, 0.01);
0168   heeEnSimRec = iB.book2D("heeEnSimRec", "", 200, 0, 0.002, 200, 0, 0.2);
0169 
0170   hefEnRec = iB.book1D("hefEnRec", "", 1000, 0, 10);
0171   hefEnSim = iB.book1D("hefEnSim", "", 1000, 0, 0.01);
0172   hefEnSimRec = iB.book2D("hefEnSimRec", "", 200, 0, 0.001, 200, 0, 0.5);
0173 
0174   hebEnRec = iB.book1D("hebEnRec", "", 1000, 0, 15);
0175   hebEnSim = iB.book1D("hebEnSim", "", 1000, 0, 0.01);
0176   hebEnSimRec = iB.book2D("hebEnSimRec", "", 200, 0, 0.02, 200, 0, 4);
0177 }
0178 
0179 void HGCalHitValidation::dqmBeginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
0180   //initiating hgc Geometry
0181   for (size_t i = 0; i < geometrySource_.size(); i++) {
0182     edm::ESHandle<HGCalDDDConstants> hgcCons = iSetup.getHandle(tok_ddd_[i]);
0183     if (hgcCons.isValid()) {
0184       hgcCons_.push_back(hgcCons.product());
0185     } else {
0186       edm::LogWarning("HGCalValid") << "Cannot initiate HGCalDDDConstants for " << geometrySource_[i] << std::endl;
0187     }
0188     edm::ESHandle<HGCalGeometry> hgcGeom = iSetup.getHandle(tok_geom_[i]);
0189     if (hgcGeom.isValid()) {
0190       hgcGeometry_.push_back(hgcGeom.product());
0191     } else {
0192       edm::LogWarning("HGCalValid") << "Cannot initiate HGCalGeometry for " << geometrySource_[i] << std::endl;
0193     }
0194   }
0195 }
0196 
0197 void HGCalHitValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0198   std::map<unsigned int, HGCHitTuple> eeHitRefs, fhHitRefs, bhHitRefs;
0199 
0200   //Accesing ee simhits
0201   edm::Handle<std::vector<PCaloHit>> eeSimHits;
0202   iEvent.getByToken(eeSimHitToken_, eeSimHits);
0203 
0204   if (eeSimHits.isValid()) {
0205     analyzeHGCalSimHit(eeSimHits, 0, heeEnSim, eeHitRefs);
0206 #ifdef EDM_ML_DEBUG
0207     for (std::map<unsigned int, HGCHitTuple>::iterator itr = eeHitRefs.begin(); itr != eeHitRefs.end(); ++itr) {
0208       int idx = std::distance(eeHitRefs.begin(), itr);
0209       edm::LogInfo("HGCalValid") << "EEHit[" << idx << "] " << std::hex << itr->first << std::dec << "; Energy "
0210                                  << std::get<0>(itr->second) << "; Position (" << std::get<1>(itr->second) << ", "
0211                                  << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ")";
0212     }
0213 #endif
0214   } else {
0215     edm::LogVerbatim("HGCalValid") << "No EE SimHit Found ";
0216   }
0217 
0218   //Accesing fh simhits
0219   edm::Handle<std::vector<PCaloHit>> fhSimHits;
0220   iEvent.getByToken(fhSimHitToken_, fhSimHits);
0221   if (fhSimHits.isValid()) {
0222     analyzeHGCalSimHit(fhSimHits, 1, hefEnSim, fhHitRefs);
0223 #ifdef EDM_ML_DEBUG
0224     for (std::map<unsigned int, HGCHitTuple>::iterator itr = fhHitRefs.begin(); itr != fhHitRefs.end(); ++itr) {
0225       int idx = std::distance(fhHitRefs.begin(), itr);
0226       edm::LogInfo("HGCalValid") << "FHHit[" << idx << "] " << std::hex << itr->first << std::dec << "; Energy "
0227                                  << std::get<0>(itr->second) << "; Position (" << std::get<1>(itr->second) << ", "
0228                                  << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ")";
0229     }
0230 #endif
0231   } else {
0232     edm::LogVerbatim("HGCalValid") << "No FH SimHit Found ";
0233   }
0234 
0235   //Accessing bh simhits
0236   edm::Handle<std::vector<PCaloHit>> bhSimHits;
0237   iEvent.getByToken(bhSimHitToken_, bhSimHits);
0238   if (bhSimHits.isValid()) {
0239     analyzeHGCalSimHit(bhSimHits, 2, hebEnSim, bhHitRefs);
0240 #ifdef EDM_ML_DEBUG
0241     for (std::map<unsigned int, HGCHitTuple>::iterator itr = bhHitRefs.begin(); itr != bhHitRefs.end(); ++itr) {
0242       int idx = std::distance(bhHitRefs.begin(), itr);
0243       edm::LogInfo("HGCalValid") << "BHHit[" << idx << "] " << std::hex << itr->first << std::dec << "; Energy "
0244                                  << std::get<0>(itr->second) << "; Position (" << std::get<1>(itr->second) << ", "
0245                                  << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ")";
0246     }
0247 #endif
0248   } else {
0249     edm::LogVerbatim("HGCalValid") << "No BH SimHit Found ";
0250   }
0251 
0252   //accessing EE Rechit information
0253   edm::Handle<HGCeeRecHitCollection> eeRecHit;
0254   iEvent.getByToken(eeRecHitToken_, eeRecHit);
0255   if (eeRecHit.isValid()) {
0256     const HGCeeRecHitCollection* theHits = (eeRecHit.product());
0257     for (auto it = theHits->begin(); it != theHits->end(); ++it) {
0258       double energy = it->energy();
0259       heeEnRec->Fill(energy);
0260       std::map<unsigned int, HGCHitTuple>::const_iterator itr = eeHitRefs.find(it->id().rawId());
0261       if (itr != eeHitRefs.end()) {
0262         GlobalPoint xyz = hgcGeometry_[0]->getPosition(it->id());
0263         heeRecVsSimX->Fill(std::get<1>(itr->second), xyz.x());
0264         heeRecVsSimY->Fill(std::get<2>(itr->second), xyz.y());
0265         heeRecVsSimZ->Fill(std::get<3>(itr->second), xyz.z());
0266         heedxVsX->Fill(std::get<1>(itr->second), (xyz.x() - std::get<1>(itr->second)));
0267         heedyVsY->Fill(std::get<2>(itr->second), (xyz.y() - std::get<2>(itr->second)));
0268         heedzVsZ->Fill(std::get<3>(itr->second), (xyz.z() - std::get<3>(itr->second)));
0269         heeEnSimRec->Fill(std::get<0>(itr->second), energy);
0270 #ifdef EDM_ML_DEBUG
0271         edm::LogInfo("HGCalValid") << "EEHit: " << std::hex << it->id().rawId() << std::dec << " Sim ("
0272                                    << std::get<0>(itr->second) << ", " << std::get<1>(itr->second) << ", "
0273                                    << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ") Rec ("
0274                                    << energy << ", " << xyz.x() << ", " << xyz.y() << ", " << xyz.z() << ")";
0275 #endif
0276       }
0277     }
0278   } else {
0279     edm::LogVerbatim("HGCalValid") << "No EE RecHit Found ";
0280   }
0281 
0282   //accessing FH Rechit information
0283   edm::Handle<HGChefRecHitCollection> fhRecHit;
0284   iEvent.getByToken(fhRecHitToken_, fhRecHit);
0285   if (fhRecHit.isValid()) {
0286     const HGChefRecHitCollection* theHits = (fhRecHit.product());
0287     for (auto it = theHits->begin(); it != theHits->end(); ++it) {
0288       double energy = it->energy();
0289       hefEnRec->Fill(energy);
0290       std::map<unsigned int, HGCHitTuple>::const_iterator itr = fhHitRefs.find(it->id().rawId());
0291       if (itr != fhHitRefs.end()) {
0292         GlobalPoint xyz = hgcGeometry_[1]->getPosition(it->id());
0293 
0294         hefRecVsSimX->Fill(std::get<1>(itr->second), xyz.x());
0295         hefRecVsSimY->Fill(std::get<2>(itr->second), xyz.y());
0296         hefRecVsSimZ->Fill(std::get<3>(itr->second), xyz.z());
0297         hefdxVsX->Fill(std::get<1>(itr->second), (xyz.x() - std::get<1>(itr->second)));
0298         hefdyVsY->Fill(std::get<2>(itr->second), (xyz.y() - std::get<2>(itr->second)));
0299         hefdzVsZ->Fill(std::get<3>(itr->second), (xyz.z() - std::get<3>(itr->second)));
0300         hefEnSimRec->Fill(std::get<0>(itr->second), energy);
0301 #ifdef EDM_ML_DEBUG
0302         edm::LogInfo("HGCalValid") << "FHHit: " << std::hex << it->id().rawId() << std::dec << " Sim ("
0303                                    << std::get<0>(itr->second) << ", " << std::get<1>(itr->second) << ", "
0304                                    << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ") Rec ("
0305                                    << energy << "," << xyz.x() << ", " << xyz.y() << ", " << xyz.z() << ")";
0306 #endif
0307       }
0308     }
0309   } else {
0310     edm::LogVerbatim("HGCalValid") << "No FH RecHit Found ";
0311   }
0312 
0313   //accessing BH Rechit information
0314   edm::Handle<HGChebRecHitCollection> bhRecHit;
0315   iEvent.getByToken(bhRecHitToken_, bhRecHit);
0316   if (bhRecHit.isValid()) {
0317     const HGChebRecHitCollection* theHits = (bhRecHit.product());
0318     for (auto it = theHits->begin(); it != theHits->end(); ++it) {
0319       double energy = it->energy();
0320       hebEnRec->Fill(energy);
0321       std::map<unsigned int, HGCHitTuple>::const_iterator itr = bhHitRefs.find(it->id().rawId());
0322       GlobalPoint xyz = hgcGeometry_[2]->getPosition(it->id());
0323       if (itr != bhHitRefs.end()) {
0324         float ang3 = xyz.phi().value();  // returns the phi in radians
0325         double fac = sinh(std::get<1>(itr->second));
0326         double pT = std::get<3>(itr->second) / fac;
0327         double xp = pT * cos(std::get<2>(itr->second));
0328         double yp = pT * sin(std::get<2>(itr->second));
0329         hebRecVsSimX->Fill(xp, xyz.x());
0330         hebRecVsSimY->Fill(yp, xyz.y());
0331         hebRecVsSimZ->Fill(std::get<3>(itr->second), xyz.z());
0332         hebdEtaVsEta->Fill(std::get<1>(itr->second), (xyz.eta() - std::get<1>(itr->second)));
0333         hebdPhiVsPhi->Fill(std::get<2>(itr->second), (ang3 - std::get<2>(itr->second)));
0334         hebdzVsZ->Fill(std::get<3>(itr->second), (xyz.z() - std::get<3>(itr->second)));
0335         hebEnSimRec->Fill(std::get<0>(itr->second), energy);
0336 
0337 #ifdef EDM_ML_DEBUG
0338         edm::LogInfo("HGCalValid") << "BHHit: " << std::hex << it->id().rawId() << std::dec << " Sim ("
0339                                    << std::get<0>(itr->second) << ", " << std::get<1>(itr->second) << ", "
0340                                    << std::get<2>(itr->second) << ", " << std::get<3>(itr->second) << ") Rec ("
0341                                    << energy << ", " << xyz.x() << ", " << xyz.y() << ", " << xyz.z() << ")\n";
0342 #endif
0343       }
0344     }
0345   } else {
0346     edm::LogVerbatim("HGCalValid") << "No BH RecHit Found ";
0347   }
0348 }
0349 
0350 void HGCalHitValidation::analyzeHGCalSimHit(edm::Handle<std::vector<PCaloHit>> const& simHits,
0351                                             int idet,
0352                                             MonitorElement* hist,
0353                                             std::map<unsigned int, HGCHitTuple>& hitRefs) {
0354   const HGCalTopology& hTopo = hgcGeometry_[idet]->topology();
0355   for (std::vector<PCaloHit>::const_iterator simHit = simHits->begin(); simHit != simHits->end(); ++simHit) {
0356     int subdet, zside, layer, wafer, celltype, cell;
0357     HGCalTestNumbering::unpackHexagonIndex(simHit->id(), subdet, zside, layer, wafer, celltype, cell);
0358     std::pair<float, float> xy = hgcCons_[idet]->locateCell(cell, layer, wafer, false);
0359     float zp = hgcCons_[idet]->waferZ(layer, false);
0360     if (zside < 0)
0361       zp = -zp;
0362     float xp = (zp < 0) ? -xy.first / 10 : xy.first / 10;
0363     float yp = xy.second / 10.0;
0364 
0365     //skip this hit if after ganging it is not valid
0366     std::pair<int, int> recoLayerCell = hgcCons_[idet]->simToReco(cell, layer, wafer, hTopo.detectorType());
0367     cell = recoLayerCell.first;
0368     layer = recoLayerCell.second;
0369 
0370     //skip this hit if after ganging it is not valid
0371     if (layer < 0 || cell < 0) {
0372     } else {
0373       //assign the RECO DetId
0374       HGCalDetId id = HGCalDetId((ForwardSubdetector)(subdet), zside, layer, celltype, wafer, cell);
0375       float energy = simHit->energy();
0376 
0377       float energySum(energy);
0378       if (hitRefs.count(id.rawId()) != 0)
0379         energySum += std::get<0>(hitRefs[id.rawId()]);
0380       hitRefs[id.rawId()] = std::make_tuple(energySum, xp, yp, zp);
0381       hist->Fill(energy);
0382     }
0383   }
0384 }
0385 
0386 //define this as a plug-in
0387 #include "FWCore/Framework/interface/MakerMacros.h"
0388 DEFINE_FWK_MODULE(HGCalHitValidation);