File indexing completed on 2025-01-22 07:34:36
0001 #include <memory>
0002
0003 #include "Validation/MuonCSCDigis/interface/CSCDigiMatcher.h"
0004
0005 using namespace std;
0006
0007 CSCDigiMatcher::CSCDigiMatcher(const edm::ParameterSet& pset, edm::ConsumesCollector&& iC) {
0008 const auto& wireDigi = pset.getParameterSet("cscWireDigi");
0009 verboseWG_ = wireDigi.getParameter<int>("verbose");
0010 minBXWire_ = wireDigi.getParameter<int>("minBX");
0011 maxBXWire_ = wireDigi.getParameter<int>("maxBX");
0012 matchDeltaWG_ = wireDigi.getParameter<int>("matchDeltaWG");
0013
0014 const auto& comparatorDigi = pset.getParameterSet("cscComparatorDigi");
0015 verboseComparator_ = comparatorDigi.getParameter<int>("verbose");
0016 minBXComparator_ = comparatorDigi.getParameter<int>("minBX");
0017 maxBXComparator_ = comparatorDigi.getParameter<int>("maxBX");
0018 matchDeltaComparator_ = comparatorDigi.getParameter<int>("matchDeltaStrip");
0019
0020 const auto& stripDigi = pset.getParameterSet("cscStripDigi");
0021 verboseStrip_ = stripDigi.getParameter<int>("verbose");
0022 minBXStrip_ = stripDigi.getParameter<int>("minBX");
0023 maxBXStrip_ = stripDigi.getParameter<int>("maxBX");
0024 matchDeltaStrip_ = stripDigi.getParameter<int>("matchDeltaStrip");
0025
0026
0027 muonSimHitMatcher_ = std::make_shared<CSCSimHitMatcher>(pset, std::move(iC));
0028
0029 comparatorDigiInput_ =
0030 iC.consumes<CSCComparatorDigiCollection>(comparatorDigi.getParameter<edm::InputTag>("inputTag"));
0031 stripDigiInput_ = iC.consumes<CSCStripDigiCollection>(stripDigi.getParameter<edm::InputTag>("inputTag"));
0032 wireDigiInput_ = iC.consumes<CSCWireDigiCollection>(wireDigi.getParameter<edm::InputTag>("inputTag"));
0033 }
0034
0035 void CSCDigiMatcher::init(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0036 muonSimHitMatcher_->init(iEvent, iSetup);
0037
0038 iEvent.getByToken(comparatorDigiInput_, comparatorDigisH_);
0039 iEvent.getByToken(stripDigiInput_, stripDigisH_);
0040 iEvent.getByToken(wireDigiInput_, wireDigisH_);
0041 }
0042
0043
0044 void CSCDigiMatcher::match(const SimTrack& t, const SimVertex& v) {
0045
0046 muonSimHitMatcher_->match(t, v);
0047
0048
0049 const CSCComparatorDigiCollection& comparators = *comparatorDigisH_.product();
0050 const CSCStripDigiCollection& strips = *stripDigisH_.product();
0051 const CSCWireDigiCollection& wires = *wireDigisH_.product();
0052
0053 clear();
0054
0055
0056 matchComparatorsToSimTrack(comparators);
0057 matchStripsToSimTrack(strips);
0058 matchWiresToSimTrack(wires);
0059 }
0060
0061 void CSCDigiMatcher::matchComparatorsToSimTrack(const CSCComparatorDigiCollection& comparators) {
0062 const auto& det_ids = muonSimHitMatcher_->detIds(0);
0063 for (const auto& id : det_ids) {
0064 CSCDetId layer_id(id);
0065
0066 const auto& hit_comparators = muonSimHitMatcher_->hitStripsInDetId(id, matchDeltaStrip_);
0067 if (verboseComparator_) {
0068 cout << "hit_comparators_fat, CSCid " << layer_id << " ";
0069 copy(hit_comparators.begin(), hit_comparators.end(), ostream_iterator<int>(cout, " "));
0070 cout << endl;
0071 }
0072
0073 int ndigis = 0;
0074 const auto& comp_digis_in_det = comparators.get(layer_id);
0075 for (auto c = comp_digis_in_det.first; c != comp_digis_in_det.second; ++c) {
0076 if (verboseComparator_)
0077 edm::LogInfo("CSCDigiMatcher") << "sdigi " << layer_id << " (comparator, comparator, Tbin ) " << *c;
0078
0079
0080 if (c->getTimeBin() < minBXComparator_ || c->getTimeBin() > maxBXComparator_)
0081 continue;
0082
0083 ndigis++;
0084
0085 int comparator = c->getStrip();
0086
0087 if (hit_comparators.find(comparator) == hit_comparators.end())
0088 continue;
0089
0090 if (verboseComparator_)
0091 edm::LogInfo("CSCDigiMatcher") << "Matched comparator " << *c;
0092 detid_to_comparators_[id].push_back(*c);
0093 chamber_to_comparators_[layer_id.chamberId().rawId()].push_back(*c);
0094 }
0095 detid_to_totalcomparators_[id] = ndigis;
0096 }
0097 }
0098
0099 void CSCDigiMatcher::matchStripsToSimTrack(const CSCStripDigiCollection& strips) {
0100 for (auto detUnitIt = strips.begin(); detUnitIt != strips.end(); ++detUnitIt) {
0101 const CSCDetId& id = (*detUnitIt).first;
0102 const auto& range = (*detUnitIt).second;
0103 for (auto digiIt = range.first; digiIt != range.second; ++digiIt) {
0104 if (id.station() == 1 and (id.ring() == 1 or id.ring() == 4))
0105 if (verboseStrip_)
0106 edm::LogInfo("CSCDigiMatcher") << "CSCid " << id << " Strip digi (strip, strip, Tbin ) " << (*digiIt);
0107 }
0108 }
0109
0110 const auto& det_ids = muonSimHitMatcher_->detIds(0);
0111 for (const auto& id : det_ids) {
0112 CSCDetId layer_id(id);
0113
0114 const auto& hit_strips = muonSimHitMatcher_->hitStripsInDetId(id, matchDeltaStrip_);
0115 if (verboseStrip_) {
0116 cout << "hit_strips_fat, CSCid " << layer_id << " ";
0117 copy(hit_strips.begin(), hit_strips.end(), ostream_iterator<int>(cout, " "));
0118 cout << endl;
0119 }
0120
0121 int ndigis = 0;
0122
0123 const auto& strip_digis_in_det = strips.get(layer_id);
0124 for (auto c = strip_digis_in_det.first; c != strip_digis_in_det.second; ++c) {
0125
0126
0127 if (verboseStrip_)
0128 edm::LogInfo("CSCDigiMatcher") << "sdigi " << layer_id << " (strip, ADC ) " << *c;
0129
0130 ndigis++;
0131
0132 int strip = c->getStrip();
0133
0134 if (hit_strips.find(strip) == hit_strips.end())
0135 continue;
0136
0137 if (verboseStrip_)
0138 edm::LogInfo("CSCDigiMatcher") << "Matched strip " << *c;
0139 detid_to_strips_[id].push_back(*c);
0140 chamber_to_strips_[layer_id.chamberId().rawId()].push_back(*c);
0141 }
0142 detid_to_totalstrips_[id] = ndigis;
0143 }
0144 }
0145
0146 void CSCDigiMatcher::matchWiresToSimTrack(const CSCWireDigiCollection& wires) {
0147 const auto& det_ids = muonSimHitMatcher_->detIds(0);
0148 for (const auto& id : det_ids) {
0149 CSCDetId layer_id(id);
0150
0151 const auto& hit_wires = muonSimHitMatcher_->hitWiregroupsInDetId(id, matchDeltaWG_);
0152 if (verboseWG_) {
0153 cout << "hit_wires ";
0154 copy(hit_wires.begin(), hit_wires.end(), ostream_iterator<int>(cout, " "));
0155 cout << endl;
0156 }
0157
0158 int ndigis = 0;
0159
0160 const auto& wire_digis_in_det = wires.get(layer_id);
0161 for (auto w = wire_digis_in_det.first; w != wire_digis_in_det.second; ++w) {
0162 if (verboseStrip_)
0163 edm::LogInfo("CSCDigiMatcher") << "wdigi " << layer_id << " (wire, Tbin ) " << *w;
0164
0165
0166 if (w->getTimeBin() < minBXWire_ || w->getTimeBin() > maxBXWire_)
0167 continue;
0168
0169 ndigis++;
0170
0171 int wg = w->getWireGroup();
0172
0173 if (hit_wires.find(wg) == hit_wires.end())
0174 continue;
0175
0176 if (verboseStrip_)
0177 edm::LogInfo("CSCDigiMatcher") << "Matched wire digi " << *w << endl;
0178 detid_to_wires_[id].push_back(*w);
0179 chamber_to_wires_[layer_id.chamberId().rawId()].push_back(*w);
0180 }
0181 detid_to_totalwires_[id] = ndigis;
0182 }
0183 }
0184
0185 std::set<unsigned int> CSCDigiMatcher::detIdsComparator(int csc_type) const {
0186 return selectDetIds(detid_to_comparators_, csc_type);
0187 }
0188
0189 std::set<unsigned int> CSCDigiMatcher::detIdsStrip(int csc_type) const {
0190 return selectDetIds(detid_to_strips_, csc_type);
0191 }
0192
0193 std::set<unsigned int> CSCDigiMatcher::detIdsWire(int csc_type) const {
0194 return selectDetIds(detid_to_wires_, csc_type);
0195 }
0196
0197 std::set<unsigned int> CSCDigiMatcher::chamberIdsComparator(int csc_type) const {
0198 return selectDetIds(chamber_to_comparators_, csc_type);
0199 }
0200
0201 std::set<unsigned int> CSCDigiMatcher::chamberIdsStrip(int csc_type) const {
0202 return selectDetIds(chamber_to_strips_, csc_type);
0203 }
0204
0205 std::set<unsigned int> CSCDigiMatcher::chamberIdsWire(int csc_type) const {
0206 return selectDetIds(chamber_to_wires_, csc_type);
0207 }
0208
0209 const CSCComparatorDigiContainer& CSCDigiMatcher::comparatorDigisInDetId(unsigned int detid) const {
0210 if (detid_to_comparators_.find(detid) == detid_to_comparators_.end())
0211 return no_comparators_;
0212 return detid_to_comparators_.at(detid);
0213 }
0214
0215 const CSCComparatorDigiContainer& CSCDigiMatcher::comparatorDigisInChamber(unsigned int detid) const {
0216 if (chamber_to_comparators_.find(detid) == chamber_to_comparators_.end())
0217 return no_comparators_;
0218 return chamber_to_comparators_.at(detid);
0219 }
0220
0221 const CSCStripDigiContainer& CSCDigiMatcher::stripDigisInDetId(unsigned int detid) const {
0222 if (detid_to_strips_.find(detid) == detid_to_strips_.end())
0223 return no_strips_;
0224 return detid_to_strips_.at(detid);
0225 }
0226
0227 const CSCStripDigiContainer& CSCDigiMatcher::stripDigisInChamber(unsigned int detid) const {
0228 if (chamber_to_strips_.find(detid) == chamber_to_strips_.end())
0229 return no_strips_;
0230 return chamber_to_strips_.at(detid);
0231 }
0232
0233 const CSCWireDigiContainer& CSCDigiMatcher::wireDigisInDetId(unsigned int detid) const {
0234 if (detid_to_wires_.find(detid) == detid_to_wires_.end())
0235 return no_wires_;
0236 return detid_to_wires_.at(detid);
0237 }
0238
0239 const CSCWireDigiContainer& CSCDigiMatcher::wireDigisInChamber(unsigned int detid) const {
0240 if (chamber_to_wires_.find(detid) == chamber_to_wires_.end())
0241 return no_wires_;
0242 return chamber_to_wires_.at(detid);
0243 }
0244
0245 int CSCDigiMatcher::nLayersWithComparatorInChamber(unsigned int detid) const {
0246 int nLayers = 0;
0247 CSCDetId chamberId(detid);
0248 for (int i = 1; i <= 6; ++i) {
0249 CSCDetId layerId(chamberId.endcap(), chamberId.station(), chamberId.ring(), chamberId.chamber(), i);
0250 if (!comparatorDigisInDetId(layerId.rawId()).empty()) {
0251 nLayers++;
0252 }
0253 }
0254 return nLayers;
0255 }
0256
0257 int CSCDigiMatcher::nLayersWithStripInChamber(unsigned int detid) const {
0258 int nLayers = 0;
0259 CSCDetId chamberId(detid);
0260 for (int i = 1; i <= 6; ++i) {
0261 CSCDetId layerId(chamberId.endcap(), chamberId.station(), chamberId.ring(), chamberId.chamber(), i);
0262 if (!stripDigisInDetId(layerId.rawId()).empty()) {
0263 nLayers++;
0264 }
0265 }
0266 return nLayers;
0267 }
0268
0269 int CSCDigiMatcher::nLayersWithWireInChamber(unsigned int detid) const {
0270 int nLayers = 0;
0271 CSCDetId chamberId(detid);
0272 for (int i = 1; i <= 6; ++i) {
0273 CSCDetId layerId(chamberId.endcap(), chamberId.station(), chamberId.ring(), chamberId.chamber(), i);
0274 if (!wireDigisInDetId(layerId.rawId()).empty()) {
0275 nLayers++;
0276 }
0277 }
0278 return nLayers;
0279 }
0280
0281 int CSCDigiMatcher::nCoincidenceComparatorChambers(int min_n_layers) const {
0282 int result = 0;
0283 const auto& chamber_ids = chamberIdsComparator();
0284 for (const auto& id : chamber_ids) {
0285 if (nLayersWithComparatorInChamber(id) >= min_n_layers)
0286 result += 1;
0287 }
0288 return result;
0289 }
0290
0291 int CSCDigiMatcher::nCoincidenceStripChambers(int min_n_layers) const {
0292 int result = 0;
0293 const auto& chamber_ids = chamberIdsStrip();
0294 for (const auto& id : chamber_ids) {
0295 if (nLayersWithStripInChamber(id) >= min_n_layers)
0296 result += 1;
0297 }
0298 return result;
0299 }
0300
0301 int CSCDigiMatcher::nCoincidenceWireChambers(int min_n_layers) const {
0302 int result = 0;
0303 const auto& chamber_ids = chamberIdsWire();
0304 for (const auto& id : chamber_ids) {
0305 if (nLayersWithWireInChamber(id) >= min_n_layers)
0306 result += 1;
0307 }
0308 return result;
0309 }
0310
0311 std::set<int> CSCDigiMatcher::comparatorsInDetId(unsigned int detid) const {
0312 set<int> result;
0313 const auto& digis = comparatorDigisInDetId(detid);
0314 for (const auto& d : digis) {
0315 result.insert(d.getHalfStrip());
0316 }
0317 return result;
0318 }
0319
0320 std::set<int> CSCDigiMatcher::stripsInDetId(unsigned int detid) const {
0321 set<int> result;
0322 const auto& digis = stripDigisInDetId(detid);
0323 for (const auto& d : digis) {
0324 result.insert(d.getStrip());
0325 }
0326 return result;
0327 }
0328
0329 std::set<int> CSCDigiMatcher::wiregroupsInDetId(unsigned int detid) const {
0330 set<int> result;
0331 const auto& digis = wireDigisInDetId(detid);
0332 for (const auto& d : digis) {
0333 result.insert(d.getWireGroup());
0334 }
0335 return result;
0336 }
0337
0338 std::set<int> CSCDigiMatcher::comparatorsInChamber(unsigned int detid, int max_gap_to_fill) const {
0339 set<int> result;
0340 const auto& digis = comparatorDigisInChamber(detid);
0341 for (const auto& d : digis) {
0342 result.insert(d.getStrip());
0343 }
0344 if (max_gap_to_fill > 0) {
0345 int prev = -111;
0346 for (const auto& s : result) {
0347 if (s - prev > 1 && s - prev - 1 <= max_gap_to_fill) {
0348 for (int fill_s = prev + 1; fill_s < s; ++fill_s)
0349 result.insert(fill_s);
0350 }
0351 prev = s;
0352 }
0353 }
0354
0355 return result;
0356 }
0357
0358 std::set<int> CSCDigiMatcher::stripsInChamber(unsigned int detid, int max_gap_to_fill) const {
0359 set<int> result;
0360 const auto& digis = stripDigisInChamber(detid);
0361 for (const auto& d : digis) {
0362 result.insert(d.getStrip());
0363 }
0364 if (max_gap_to_fill > 0) {
0365 int prev = -111;
0366 for (const auto& s : result) {
0367 if (s - prev > 1 && s - prev - 1 <= max_gap_to_fill) {
0368 for (int fill_s = prev + 1; fill_s < s; ++fill_s)
0369 result.insert(fill_s);
0370 }
0371 prev = s;
0372 }
0373 }
0374
0375 return result;
0376 }
0377
0378 std::set<int> CSCDigiMatcher::wiregroupsInChamber(unsigned int detid, int max_gap_to_fill) const {
0379 set<int> result;
0380 const auto& digis = wireDigisInChamber(detid);
0381 for (const auto& d : digis) {
0382 result.insert(d.getWireGroup());
0383 }
0384 if (max_gap_to_fill > 0) {
0385 int prev = -111;
0386 for (const auto& w : result) {
0387 if (w - prev > 1 && w - prev - 1 <= max_gap_to_fill) {
0388 for (int fill_w = prev + 1; fill_w < w; ++fill_w)
0389 result.insert(fill_w);
0390 }
0391 prev = w;
0392 }
0393 }
0394 return result;
0395 }
0396
0397 int CSCDigiMatcher::totalComparators(unsigned int detid) const {
0398 if (detid_to_totalcomparators_.find(detid) == detid_to_totalcomparators_.end())
0399 return 0;
0400 return detid_to_totalcomparators_.at(detid);
0401 }
0402
0403 int CSCDigiMatcher::totalStrips(unsigned int detid) const {
0404 if (detid_to_totalstrips_.find(detid) == detid_to_totalstrips_.end())
0405 return 0;
0406 return detid_to_totalstrips_.at(detid);
0407 }
0408
0409 int CSCDigiMatcher::totalWires(unsigned int detid) const {
0410 if (detid_to_totalwires_.find(detid) == detid_to_totalwires_.end())
0411 return 0;
0412 return detid_to_totalwires_.at(detid);
0413 }
0414
0415 void CSCDigiMatcher::clear() {
0416 detid_to_comparators_.clear();
0417 chamber_to_comparators_.clear();
0418
0419 detid_to_strips_.clear();
0420 chamber_to_strips_.clear();
0421
0422 detid_to_wires_.clear();
0423 chamber_to_wires_.clear();
0424 }