File indexing completed on 2024-09-11 04:33:26
0001 #include "Validation/MuonGEMDigis/plugins/GEMPadDigiValidation.h"
0002
0003 GEMPadDigiValidation::GEMPadDigiValidation(const edm::ParameterSet& pset)
0004 : GEMBaseValidation(pset, "GEMPadDigiValidation") {
0005 const auto& pad_pset = pset.getParameterSet("gemPadDigi");
0006 const auto& pad_tag = pad_pset.getParameter<edm::InputTag>("inputTag");
0007 pad_token_ = consumes<GEMPadDigiCollection>(pad_tag);
0008
0009 const auto& simhit_pset = pset.getParameterSet("gemSimHit");
0010 const auto& simhit_tag = simhit_pset.getParameter<edm::InputTag>("inputTag");
0011 simhit_token_ = consumes<edm::PSimHitContainer>(simhit_tag);
0012
0013 const auto& digisimlink_tag = pset.getParameter<edm::InputTag>("gemDigiSimLink");
0014 digisimlink_token_ = consumes<edm::DetSetVector<GEMDigiSimLink>>(digisimlink_tag);
0015
0016 geomToken_ = esConsumes<GEMGeometry, MuonGeometryRecord>();
0017 geomTokenBeginRun_ = esConsumes<GEMGeometry, MuonGeometryRecord, edm::Transition::BeginRun>();
0018 }
0019
0020 void GEMPadDigiValidation::bookHistograms(DQMStore::IBooker& booker,
0021 edm::Run const& Run,
0022 edm::EventSetup const& setup) {
0023 const auto gemH = setup.getHandle(geomTokenBeginRun_);
0024 if (!gemH.isValid()) {
0025 edm::LogError(kLogCategory_) << "Failed to initialize GEM geometry.";
0026 return;
0027 }
0028 const GEMGeometry* gem = gemH.product();
0029
0030
0031 booker.setCurrentFolder("GEM/Pad");
0032
0033 for (const auto& region : gem->regions()) {
0034 Int_t region_id = region->region();
0035
0036 if (detail_plot_)
0037 me_detail_occ_zr_[region_id] = bookZROccupancy(booker, region_id, "pad", "Pad");
0038
0039 for (const auto& station : region->stations()) {
0040 Int_t station_id = station->station();
0041 ME2IdsKey key2(region_id, station_id);
0042
0043 if (detail_plot_) {
0044 me_detail_occ_det_[key2] = bookDetectorOccupancy(booker, key2, station, "pad", "Pad");
0045 me_detail_pad_occ_det_[key2] = bookDetectorOccupancy(booker, key2, station, "sim_matched", "Matched Pad");
0046 }
0047
0048 const auto& superChamberVec = station->superChambers();
0049 if (superChamberVec.empty()) {
0050 edm::LogError(kLogCategory_) << "Super chambers missing for region = " << region_id
0051 << " and station = " << station_id;
0052 continue;
0053 }
0054 const GEMSuperChamber* super_chamber = superChamberVec.front();
0055 if (super_chamber == nullptr) {
0056 edm::LogError(kLogCategory_) << "Failed to find super chamber for region = " << region_id
0057 << " and station = " << station_id;
0058 continue;
0059 }
0060 for (const auto& chamber : super_chamber->chambers()) {
0061 Int_t layer_id = chamber->id().layer();
0062 ME3IdsKey key3(region_id, station_id, layer_id);
0063
0064 const auto& etaPartitionVec = chamber->etaPartitions();
0065 if (etaPartitionVec.empty() || etaPartitionVec.front() == nullptr) {
0066 edm::LogError(kLogCategory_)
0067 << "Eta partition missing or null for region, station, super chamber, chamber = (" << region_id << ", "
0068 << station_id << ", " << super_chamber->id() << ", " << chamber->id() << ")";
0069 continue;
0070 }
0071 Int_t num_pads = etaPartitionVec.front()->npads();
0072
0073 me_occ_total_pad_[key3] =
0074 bookHist1D(booker, key3, "total_pads_per_event", "Number of pad digis per event", 50, 0, 50);
0075
0076 me_pad_occ_eta_[key3] = bookHist1D(booker,
0077 key3,
0078 "sim_matched_occ_eta",
0079 "Matched Pad Eta Occupancy",
0080 16,
0081 eta_range_[station_id * 2 + 0],
0082 eta_range_[station_id * 2 + 1],
0083 "#eta");
0084
0085 me_pad_occ_phi_[key3] =
0086 bookHist1D(booker, key3, "sim_matched_occ_phi", "Matched Pad Phi Occupancy", 36, -5, 355, "#phi [degrees]");
0087
0088 if (detail_plot_) {
0089 me_detail_occ_xy_[key3] = bookXYOccupancy(booker, key3, "pad", "Pad");
0090
0091 me_detail_occ_phi_pad_[key3] = bookHist2D(booker,
0092 key3,
0093 "occ_phi_pad",
0094 "Pad Occupancy",
0095 280,
0096 -M_PI,
0097 M_PI,
0098 num_pads / 2,
0099 0,
0100 num_pads,
0101 "#phi [rad]",
0102 "Pad number");
0103
0104 me_detail_occ_pad_[key3] =
0105 bookHist1D(booker, key3, "occ_pad", "Pad Occupancy", num_pads, 0, num_pads, "Pad number");
0106 }
0107 }
0108 }
0109 }
0110
0111
0112 if (detail_plot_) {
0113 for (const auto& region : gem->regions()) {
0114 Int_t region_id = region->region();
0115 for (const auto& station : region->stations()) {
0116 Int_t station_id = station->station();
0117
0118 const auto& superChamberVec = station->superChambers();
0119 if (superChamberVec.empty()) {
0120 edm::LogError(kLogCategory_) << "Super chambers missing for region = " << region_id
0121 << " and station = " << station_id;
0122 continue;
0123 }
0124 const GEMSuperChamber* super_chamber = superChamberVec.front();
0125 if (super_chamber == nullptr) {
0126 edm::LogError(kLogCategory_) << "Failed to find super chamber for region = " << region_id
0127 << " and station = " << station_id;
0128 continue;
0129 }
0130 for (const auto& chamber : super_chamber->chambers()) {
0131 Int_t layer_id = chamber->id().layer();
0132 ME3IdsKey key3(region_id, station_id, layer_id);
0133
0134 me_detail_bx_[key3] = bookHist1D(booker, key3, "bx", "Pad Bunch Crossing", 5, -2, 3, "Bunch crossing");
0135 }
0136 }
0137 }
0138 }
0139 }
0140
0141 GEMPadDigiValidation::~GEMPadDigiValidation() {}
0142
0143 void GEMPadDigiValidation::analyze(const edm::Event& event, const edm::EventSetup& setup) {
0144 const auto& gemH = setup.getHandle(geomToken_);
0145 if (!gemH.isValid()) {
0146 edm::LogError(kLogCategory_) << "Failed to initialize GEM geometry.";
0147 return;
0148 }
0149 const GEMGeometry* gem = gemH.product();
0150
0151 edm::Handle<GEMPadDigiCollection> collection;
0152 event.getByToken(pad_token_, collection);
0153 if (not collection.isValid()) {
0154 edm::LogError(kLogCategory_) << "Cannot get pads by label GEMPadToken.";
0155 return;
0156 }
0157
0158 edm::Handle<edm::DetSetVector<GEMDigiSimLink>> digiSimLink;
0159 event.getByToken(digisimlink_token_, digiSimLink);
0160 if (not digiSimLink.isValid()) {
0161 edm::LogError(kLogCategory_) << "Failed to get GEMDigiSimLink." << std::endl;
0162 return;
0163 }
0164
0165 edm::Handle<edm::PSimHitContainer> simhit_container;
0166 event.getByToken(simhit_token_, simhit_container);
0167 if (not simhit_container.isValid()) {
0168 edm::LogError(kLogCategory_) << "Failed to get PSimHitContainer." << std::endl;
0169 return;
0170 }
0171
0172 std::map<ME3IdsKey, Int_t> total_pad;
0173 for (const auto& pad_pair : *collection) {
0174 GEMDetId gemid = pad_pair.first;
0175 const auto& range = pad_pair.second;
0176
0177 if (gem->idToDet(gemid) == nullptr) {
0178 edm::LogError(kLogCategory_) << "Getting DetId failed. Discard this gem pad hit. "
0179 << "Maybe it comes from unmatched geometry." << std::endl;
0180 continue;
0181 }
0182
0183 const GEMEtaPartition* roll = gem->etaPartition(gemid);
0184 const BoundPlane& surface = roll->surface();
0185
0186 Int_t region_id = gemid.region();
0187 Int_t station_id = gemid.station();
0188 Int_t layer_id = gemid.layer();
0189 Int_t chamber_id = gemid.chamber();
0190 Int_t ieta = gemid.ieta();
0191 Int_t num_layers = gemid.nlayers();
0192
0193 ME2IdsKey key2(region_id, station_id);
0194 ME3IdsKey key3(region_id, station_id, layer_id);
0195
0196 for (auto digi = range.first; digi != range.second; ++digi) {
0197 total_pad[key3]++;
0198
0199 Int_t pad = digi->pad();
0200 Int_t bx = digi->bx();
0201
0202 const LocalPoint& local_pos = roll->centreOfPad(pad);
0203 const GlobalPoint& global_pos = surface.toGlobal(local_pos);
0204
0205 Float_t g_r = global_pos.perp();
0206 Float_t g_phi = global_pos.phi();
0207 Float_t g_x = global_pos.x();
0208 Float_t g_y = global_pos.y();
0209 Float_t g_abs_z = std::fabs(global_pos.z());
0210
0211 Int_t bin_x = getDetOccBinX(num_layers, chamber_id, layer_id);
0212
0213 if (detail_plot_) {
0214 me_detail_occ_zr_[region_id]->Fill(g_abs_z, g_r);
0215 me_detail_occ_xy_[key3]->Fill(g_x, g_y);
0216 me_detail_occ_phi_pad_[key3]->Fill(g_phi, pad);
0217 me_detail_occ_pad_[key3]->Fill(pad);
0218 me_detail_occ_det_[key2]->Fill(bin_x, ieta);
0219 me_detail_bx_[key3]->Fill(bx);
0220 }
0221 }
0222 }
0223
0224 for (const auto& region : gem->regions()) {
0225 Int_t region_id = region->region();
0226 for (const auto& station : region->stations()) {
0227 Int_t station_id = station->station();
0228 const auto& superChamberVec = station->superChambers();
0229 if (superChamberVec.empty()) {
0230 edm::LogError(kLogCategory_) << "Super chambers missing for region = " << region_id
0231 << " and station = " << station_id;
0232 continue;
0233 }
0234 const GEMSuperChamber* super_chamber = superChamberVec.front();
0235 if (super_chamber == nullptr) {
0236 edm::LogError(kLogCategory_) << "Failed to find super chamber for region = " << region_id
0237 << " and station = " << station_id;
0238 continue;
0239 }
0240 for (const auto& chamber : super_chamber->chambers()) {
0241 Int_t layer_id = chamber->id().layer();
0242 ME3IdsKey key3{region_id, station_id, layer_id};
0243 me_occ_total_pad_[key3]->Fill(total_pad[key3]);
0244 }
0245 }
0246 }
0247
0248
0249 for (const auto& simhit : *simhit_container.product()) {
0250 if (not isMuonSimHit(simhit))
0251 continue;
0252 if (gem->idToDet(simhit.detUnitId()) == nullptr) {
0253 edm::LogError(kLogCategory_) << "SimHit did not match with GEMGeometry." << std::endl;
0254 continue;
0255 }
0256
0257 GEMDetId simhit_gemid(simhit.detUnitId());
0258
0259 Int_t region_id = simhit_gemid.region();
0260 Int_t station_id = simhit_gemid.station();
0261 Int_t layer_id = simhit_gemid.layer();
0262 Int_t chamber_id = simhit_gemid.chamber();
0263 Int_t ieta = simhit_gemid.ieta();
0264 Int_t num_layers = simhit_gemid.nlayers();
0265
0266 ME2IdsKey key2{region_id, station_id};
0267 ME3IdsKey key3{region_id, station_id, layer_id};
0268
0269 const GEMEtaPartition* roll = gem->etaPartition(simhit_gemid);
0270
0271 const auto& simhit_local_pos = simhit.localPosition();
0272 const auto& simhit_global_pos = roll->surface().toGlobal(simhit_local_pos);
0273
0274 Float_t simhit_g_eta = std::abs(simhit_global_pos.eta());
0275 Float_t simhit_g_phi = toDegree(simhit_global_pos.phi());
0276
0277 auto simhit_trackId = simhit.trackId();
0278
0279 Int_t bin_x = getDetOccBinX(num_layers, chamber_id, layer_id);
0280
0281 auto links = digiSimLink->find(simhit_gemid);
0282 if (links == digiSimLink->end())
0283 continue;
0284
0285 Int_t simhit_strip = -1;
0286 for (const auto& link : *links) {
0287 if (simhit_trackId == link.getTrackId()) {
0288 simhit_strip = link.getStrip();
0289 break;
0290 }
0291 }
0292 Int_t simhit_pad = roll->padOfStrip(simhit_strip);
0293 auto range = collection->get(simhit_gemid);
0294 for (auto pad = range.first; pad != range.second; ++pad) {
0295 if (pad->pad() == simhit_pad) {
0296 me_pad_occ_eta_[key3]->Fill(simhit_g_eta);
0297 me_pad_occ_phi_[key3]->Fill(simhit_g_phi);
0298 if (detail_plot_) {
0299 me_detail_pad_occ_det_[key2]->Fill(bin_x, ieta);
0300 }
0301 break;
0302 }
0303 }
0304 }
0305 }