MtdValidation/plugins/BtlSimHitsValidation.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    Validation/MtdValidation
0004 // Class:      BtlSimHitsValidation
0005 //
0006 /**\class BtlSimHitsValidation BtlSimHitsValidation.cc Validation/MtdValidation/plugins/BtlSimHitsValidation.cc
0007
0008  Description: BTL SIM hits validation
0009
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013
0014 #include <string>
0015
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/Framework/interface/Event.h"
0018 #include "FWCore/Framework/interface/MakerMacros.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020
0021 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0022 #include "DQMServices/Core/interface/DQMStore.h"
0023
0024 #include "DataFormats/Common/interface/ValidHandle.h"
0025 #include "DataFormats/Math/interface/GeantUnits.h"
0026 #include "DataFormats/ForwardDetId/interface/BTLDetId.h"
0027
0028 #include "SimDataFormats/CrossingFrame/interface/CrossingFrame.h"
0029 #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
0030 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0031
0032 #include "Geometry/Records/interface/MTDDigiGeometryRecord.h"
0033 #include "Geometry/Records/interface/MTDTopologyRcd.h"
0034 #include "Geometry/MTDGeometryBuilder/interface/MTDGeometry.h"
0035 #include "Geometry/MTDNumberingBuilder/interface/MTDTopology.h"
0036
0037 #include "Geometry/MTDGeometryBuilder/interface/ProxyMTDTopology.h"
0038 #include "Geometry/MTDGeometryBuilder/interface/RectangularMTDTopology.h"
0039
0040 #include "Geometry/MTDCommonData/interface/MTDTopologyMode.h"
0041
0042 #include "MTDHit.h"
0043
0044 class BtlSimHitsValidation : public DQMEDAnalyzer {
0045 public:
0046   explicit BtlSimHitsValidation(const edm::ParameterSet&);
0047   ~BtlSimHitsValidation() override;
0048
0049   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0050
0051 private:
0052   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0053
0054   void analyze(const edm::Event&, const edm::EventSetup&) override;
0055
0056   // ------------ member data ------------
0057
0058   const std::string folder_;
0059   const float hitMinEnergy_;
0060
0061   edm::EDGetTokenT<CrossingFrame<PSimHit> > btlSimHitsToken_;
0062
0063   edm::ESGetToken<MTDGeometry, MTDDigiGeometryRecord> mtdgeoToken_;
0064   edm::ESGetToken<MTDTopology, MTDTopologyRcd> mtdtopoToken_;
0065
0066   // --- histograms declaration
0067
0068   MonitorElement* meNevents_;
0069
0070   MonitorElement* meNhits_;
0071   MonitorElement* meNtrkPerCell_;
0072
0073   MonitorElement* meHitEnergy_;
0074   MonitorElement* meHitLogEnergy_;
0075   MonitorElement* meHitTime_;
0076
0077   MonitorElement* meHitXlocal_;
0078   MonitorElement* meHitYlocal_;
0079   MonitorElement* meHitZlocal_;
0080
0081   MonitorElement* meOccupancy_;
0082
0083   MonitorElement* meHitX_;
0084   MonitorElement* meHitY_;
0085   MonitorElement* meHitZ_;
0086   MonitorElement* meHitPhi_;
0087   MonitorElement* meHitEta_;
0088
0089   MonitorElement* meHitTvsE_;
0090   MonitorElement* meHitEvsPhi_;
0091   MonitorElement* meHitEvsEta_;
0092   MonitorElement* meHitEvsZ_;
0093   MonitorElement* meHitTvsPhi_;
0094   MonitorElement* meHitTvsEta_;
0095   MonitorElement* meHitTvsZ_;
0096 };
0097
0098 // ------------ constructor and destructor --------------
0099 BtlSimHitsValidation::BtlSimHitsValidation(const edm::ParameterSet& iConfig)
0100     : folder_(iConfig.getParameter<std::string>("folder")),
0101       hitMinEnergy_(iConfig.getParameter<double>("hitMinimumEnergy")) {
0102   btlSimHitsToken_ = consumes<CrossingFrame<PSimHit> >(iConfig.getParameter<edm::InputTag>("inputTag"));
0103   mtdgeoToken_ = esConsumes<MTDGeometry, MTDDigiGeometryRecord>();
0104   mtdtopoToken_ = esConsumes<MTDTopology, MTDTopologyRcd>();
0105 }
0106
0107 BtlSimHitsValidation::~BtlSimHitsValidation() {}
0108
0109 // ------------ method called for each event  ------------
0110 void BtlSimHitsValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0111   using namespace edm;
0112   using namespace geant_units::operators;
0113
0114   auto geometryHandle = iSetup.getTransientHandle(mtdgeoToken_);
0115   const MTDGeometry* geom = geometryHandle.product();
0116
0117   auto topologyHandle = iSetup.getTransientHandle(mtdtopoToken_);
0118   const MTDTopology* topology = topologyHandle.product();
0119
0120   auto btlSimHitsHandle = makeValid(iEvent.getHandle(btlSimHitsToken_));
0121   MixCollection<PSimHit> btlSimHits(btlSimHitsHandle.product());
0122
0123   std::unordered_map<uint32_t, MTDHit> m_btlHits;
0124   std::unordered_map<uint32_t, std::set<int> > m_btlTrkPerCell;
0125
0126   // --- Loop over the BLT SIM hits
0127   for (auto const& simHit : btlSimHits) {
0128     // --- Use only hits compatible with the in-time bunch-crossing
0129     if (simHit.tof() < 0 || simHit.tof() > 25.)
0130       continue;
0131
0132     DetId id = simHit.detUnitId();
0133
0134     m_btlTrkPerCell[id.rawId()].insert(simHit.trackId());
0135
0136     auto simHitIt = m_btlHits.emplace(id.rawId(), MTDHit()).first;
0137
0138     // --- Accumulate the energy (in MeV) of SIM hits in the same detector cell
0139     (simHitIt->second).energy += convertUnitsTo(0.001_MeV, simHit.energyLoss());
0140
0141     // --- Get the time of the first SIM hit in the cell
0142     if ((simHitIt->second).time == 0 || simHit.tof() < (simHitIt->second).time) {
0143       (simHitIt->second).time = simHit.tof();
0144
0145       auto hit_pos = simHit.entryPoint();
0146       (simHitIt->second).x = hit_pos.x();
0147       (simHitIt->second).y = hit_pos.y();
0148       (simHitIt->second).z = hit_pos.z();
0149     }
0150
0151   }  // simHit loop
0152
0153   // ==============================================================================
0154   //  Histogram filling
0155   // ==============================================================================
0156
0157   if (!m_btlHits.empty())
0158     meNhits_->Fill(log10(m_btlHits.size()));
0159
0160   for (auto const& hit : m_btlTrkPerCell)
0161     meNtrkPerCell_->Fill((hit.second).size());
0162
0163   for (auto const& hit : m_btlHits) {
0164     meHitLogEnergy_->Fill(log10((hit.second).energy));
0165
0166     if ((hit.second).energy < hitMinEnergy_)
0167       continue;
0168
0169     // --- Get the SIM hit global position
0170     BTLDetId detId(hit.first);
0171     DetId geoId = detId.geographicalId(MTDTopologyMode::crysLayoutFromTopoMode(topology->getMTDTopologyMode()));
0172     const MTDGeomDet* thedet = geom->idToDet(geoId);
0173     if (thedet == nullptr)
0174       throw cms::Exception("BtlSimHitsValidation") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0175                                                    << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0176     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0177     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0178
0179     Local3DPoint local_point(
0180         convertMmToCm((hit.second).x), convertMmToCm((hit.second).y), convertMmToCm((hit.second).z));
0181
0182     local_point = topo.pixelToModuleLocalPoint(local_point, detId.row(topo.nrows()), detId.column(topo.nrows()));
0183     const auto& global_point = thedet->toGlobal(local_point);
0184
0185     // --- Fill the histograms
0186     meHitEnergy_->Fill((hit.second).energy);
0187     meHitTime_->Fill((hit.second).time);
0188
0189     meHitXlocal_->Fill((hit.second).x);
0190     meHitYlocal_->Fill((hit.second).y);
0191     meHitZlocal_->Fill((hit.second).z);
0192
0193     meOccupancy_->Fill(global_point.z(), global_point.phi());
0194
0195     meHitX_->Fill(global_point.x());
0196     meHitY_->Fill(global_point.y());
0197     meHitZ_->Fill(global_point.z());
0198     meHitPhi_->Fill(global_point.phi());
0199     meHitEta_->Fill(global_point.eta());
0200
0201     meHitTvsE_->Fill((hit.second).energy, (hit.second).time);
0202     meHitEvsPhi_->Fill(global_point.phi(), (hit.second).energy);
0203     meHitEvsEta_->Fill(global_point.eta(), (hit.second).energy);
0204     meHitEvsZ_->Fill(global_point.z(), (hit.second).energy);
0205     meHitTvsPhi_->Fill(global_point.phi(), (hit.second).time);
0206     meHitTvsEta_->Fill(global_point.eta(), (hit.second).time);
0207     meHitTvsZ_->Fill(global_point.z(), (hit.second).time);
0208
0209   }  // hit loop
0210
0211   // --- This is to count the number of processed events, needed in the harvesting step
0212   meNevents_->Fill(0.5);
0213 }
0214
0215 // ------------ method for histogram booking ------------
0216 void BtlSimHitsValidation::bookHistograms(DQMStore::IBooker& ibook,
0217                                           edm::Run const& run,
0218                                           edm::EventSetup const& iSetup) {
0219   ibook.setCurrentFolder(folder_);
0220
0221   // --- histograms booking
0222
0223   meNevents_ = ibook.book1D("BtlNevents", "Number of events", 1, 0., 1.);
0224
0225   meNhits_ = ibook.book1D("BtlNhits", "Number of BTL cells with SIM hits;log_{10}(N_{BTL cells})", 100, 0., 5.25);
0226   meNtrkPerCell_ = ibook.book1D("BtlNtrkPerCell", "Number of tracks per BTL cell;N_{trk}", 10, 0., 10.);
0227
0228   meHitEnergy_ = ibook.book1D("BtlHitEnergy", "BTL SIM hits energy;E_{SIM} [MeV]", 100, 0., 20.);
0229   meHitLogEnergy_ = ibook.book1D("BtlHitLogEnergy", "BTL SIM hits energy;log_{10}(E_{SIM} [MeV])", 200, -6., 3.);
0230   meHitTime_ = ibook.book1D("BtlHitTime", "BTL SIM hits ToA;ToA_{SIM} [ns]", 100, 0., 25.);
0231
0232   meHitXlocal_ = ibook.book1D("BtlHitXlocal", "BTL SIM local X;X_{SIM}^{LOC} [mm]", 100, -30., 30.);
0233   meHitYlocal_ = ibook.book1D("BtlHitYlocal", "BTL SIM local Y;Y_{SIM}^{LOC} [mm]", 100, -1.65, 1.65);
0234   meHitZlocal_ = ibook.book1D("BtlHitZlocal", "BTL SIM local z;z_{SIM}^{LOC} [mm]", 100, -2., 2.);
0235
0236   meOccupancy_ = ibook.book2D(
0237       "BtlOccupancy", "BTL SIM hits occupancy;z_{SIM} [cm];#phi_{SIM} [rad]", 130, -260., 260., 200, -3.15, 3.15);
0238
0239   meHitX_ = ibook.book1D("BtlHitX", "BTL SIM hits X;X_{SIM} [cm]", 100, -120., 120.);
0240   meHitY_ = ibook.book1D("BtlHitY", "BTL SIM hits Y;Y_{SIM} [cm]", 100, -120., 120.);
0241   meHitZ_ = ibook.book1D("BtlHitZ", "BTL SIM hits Z;Z_{SIM} [cm]", 100, -260., 260.);
0242   meHitPhi_ = ibook.book1D("BtlHitPhi", "BTL SIM hits #phi;#phi_{SIM} [rad]", 200, -3.15, 3.15);
0243   meHitEta_ = ibook.book1D("BtlHitEta", "BTL SIM hits #eta;#eta_{SIM}", 100, -1.55, 1.55);
0244
0245   meHitTvsE_ =
0246       ibook.bookProfile("BtlHitTvsE", "BTL SIM time vs energy;E_{SIM} [MeV];T_{SIM} [ns]", 50, 0., 20., 0., 100.);
0247   meHitEvsPhi_ = ibook.bookProfile(
0248       "BtlHitEvsPhi", "BTL SIM energy vs #phi;#phi_{SIM} [rad];E_{SIM} [MeV]", 50, -3.15, 3.15, 0., 100.);
0249   meHitEvsEta_ =
0250       ibook.bookProfile("BtlHitEvsEta", "BTL SIM energy vs #eta;#eta_{SIM};E_{SIM} [MeV]", 50, -1.55, 1.55, 0., 100.);
0251   meHitEvsZ_ =
0252       ibook.bookProfile("BtlHitEvsZ", "BTL SIM energy vs Z;Z_{SIM} [cm];E_{SIM} [MeV]", 50, -260., 260., 0., 100.);
0253   meHitTvsPhi_ = ibook.bookProfile(
0254       "BtlHitTvsPhi", "BTL SIM time vs #phi;#phi_{SIM} [rad];T_{SIM} [ns]", 50, -3.15, 3.15, 0., 100.);
0255   meHitTvsEta_ =
0256       ibook.bookProfile("BtlHitTvsEta", "BTL SIM time vs #eta;#eta_{SIM};T_{SIM} [ns]", 50, -1.55, 1.55, 0., 100.);
0257   meHitTvsZ_ =
0258       ibook.bookProfile("BtlHitTvsZ", "BTL SIM time vs Z;Z_{SIM} [cm];T_{SIM} [ns]", 50, -260., 260., 0., 100.);
0259 }
0260
0261 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0262 void BtlSimHitsValidation::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0263   edm::ParameterSetDescription desc;
0264
0265   desc.add<std::string>("folder", "MTD/BTL/SimHits");
0266   desc.add<edm::InputTag>("inputTag", edm::InputTag("mix", "g4SimHitsFastTimerHitsBarrel"));
0267   desc.add<double>("hitMinimumEnergy", 1.);  // [MeV]
0268
0269   descriptions.add("btlSimHitsValid", desc);
0270 }
0271
0272 DEFINE_FWK_MODULE(BtlSimHitsValidation);