File indexing completed on 2024-09-07 04:38:10
0001
0002
0003
0004
0005
0006
0007
0008 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/Framework/interface/ESHandle.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013
0014 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0015 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0016 #include "Geometry/GEMGeometry/interface/GEMGeometry.h"
0017 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0018 #include "DataFormats/MuonDetId/interface/GEMDetId.h"
0019 #include "DataFormats/GEMDigi/interface/GEMDigiCollection.h"
0020 #include "DataFormats/Common/interface/DetSetVector.h"
0021 #include "SimDataFormats/TrackerDigiSimLink/interface/StripDigiSimLink.h"
0022 #include <map>
0023 #include <set>
0024
0025 #include "DataFormats/Common/interface/DetSet.h"
0026
0027 #include <iostream>
0028
0029 using namespace std;
0030
0031 class GEMDigiReader : public edm::one::EDAnalyzer<> {
0032 public:
0033 explicit GEMDigiReader(const edm::ParameterSet& pset);
0034
0035 ~GEMDigiReader() override {}
0036
0037 void analyze(const edm::Event&, const edm::EventSetup&) override;
0038
0039 private:
0040 const edm::EDGetTokenT<edm::PSimHitContainer> simhitToken_;
0041 const edm::EDGetTokenT<GEMDigiCollection> gemDigiToken_;
0042 const edm::EDGetTokenT<edm::DetSetVector<StripDigiSimLink> > gemDigiSimLinkToken_;
0043 const edm::ESGetToken<GEMGeometry, MuonGeometryRecord> geomToken_;
0044 };
0045
0046 GEMDigiReader::GEMDigiReader(const edm::ParameterSet& pset)
0047 : simhitToken_(consumes<edm::PSimHitContainer>(pset.getParameter<edm::InputTag>("simhitToken"))),
0048 gemDigiToken_(consumes<GEMDigiCollection>(pset.getParameter<edm::InputTag>("gemDigiToken"))),
0049 gemDigiSimLinkToken_(
0050 consumes<edm::DetSetVector<StripDigiSimLink> >(pset.getParameter<edm::InputTag>("gemDigiSimLinkToken"))),
0051 geomToken_(esConsumes<GEMGeometry, MuonGeometryRecord>()) {}
0052
0053 void GEMDigiReader::analyze(const edm::Event& event, const edm::EventSetup& eventSetup) {
0054 LogDebug("GEMDigiReader") << "--- Run: " << event.id().run() << " Event: " << event.id().event() << endl;
0055
0056 const auto& pDD = eventSetup.getHandle(geomToken_);
0057
0058 const auto& simHits = event.getHandle(simhitToken_);
0059
0060 const auto& digis = event.getHandle(gemDigiToken_);
0061
0062 const auto& thelinkDigis = event.getHandle(gemDigiSimLinkToken_);
0063
0064 GEMDigiCollection::DigiRangeIterator detUnitIt;
0065 for (detUnitIt = digis->begin(); detUnitIt != digis->end(); ++detUnitIt) {
0066 const GEMDetId& id = (*detUnitIt).first;
0067 const GEMEtaPartition* roll = pDD->etaPartition(id);
0068
0069
0070
0071
0072 LogDebug("GEMDigiReader") << "--------------" << endl;
0073 LogDebug("GEMDigiReader") << "id: " << id.rawId() << " number of strips " << roll->nstrips() << endl;
0074
0075
0076 const GEMDigiCollection::Range& range = (*detUnitIt).second;
0077 for (GEMDigiCollection::const_iterator digiIt = range.first; digiIt != range.second; ++digiIt) {
0078 LogDebug("GEMDigiReader") << " digi " << *digiIt << endl;
0079 if (digiIt->strip() < 1 || digiIt->strip() > roll->nstrips()) {
0080 LogDebug("GEMDigiReader") << " XXXXXXXXXXXXX Problemt with " << id
0081 << " a digi has strip# = " << digiIt->strip() << endl;
0082 }
0083 for (const auto& simHit : *simHits) {
0084 GEMDetId rpcId(simHit.detUnitId());
0085 if (rpcId == id && abs(simHit.particleType()) == 13) {
0086 LogDebug("GEMDigiReader") << "entry: " << simHit.entryPoint() << endl
0087 << "exit: " << simHit.exitPoint() << endl
0088 << "TOF: " << simHit.timeOfFlight() << endl;
0089 }
0090 }
0091 }
0092 }
0093
0094 for (edm::DetSetVector<StripDigiSimLink>::const_iterator itlink = thelinkDigis->begin();
0095 itlink != thelinkDigis->end();
0096 itlink++) {
0097 for (edm::DetSet<StripDigiSimLink>::const_iterator link_iter = itlink->data.begin();
0098 link_iter != itlink->data.end();
0099 ++link_iter) {
0100 int detid = itlink->detId();
0101 int ev = link_iter->eventId().event();
0102 float frac = link_iter->fraction();
0103 int strip = link_iter->channel();
0104 int trkid = link_iter->SimTrackId();
0105 int bx = link_iter->eventId().bunchCrossing();
0106 LogDebug("GEMDigiReader") << "DetUnit: " << GEMDetId(detid) << " Event ID: " << ev << " trkId: " << trkid
0107 << " Strip: " << strip << " Bx: " << bx << " frac: " << frac << endl;
0108 }
0109 }
0110
0111 LogDebug("GEMDigiReader") << "--------------" << endl;
0112 }
0113
0114 #include "FWCore/Framework/interface/MakerMacros.h"
0115 DEFINE_FWK_MODULE(GEMDigiReader);