File indexing completed on 2024-04-06 12:20:48
0001 #include <algorithm>
0002 #include "L1Trigger/L1TMuon/interface/deprecate/DTBunchCrossingCleaner.h"
0003 #include "L1Trigger/L1TMuon/interface/MuonTriggerPrimitive.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008
0009 #include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhDigi.h"
0010 #include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambThDigi.h"
0011 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0012
0013 using namespace L1TMuon;
0014
0015 namespace {
0016 typedef edm::ParameterSet PSet;
0017 }
0018
0019 DTBunchCrossingCleaner::DTBunchCrossingCleaner(const PSet& ps) : bx_window_size(ps.getParameter<int>("bxWindowSize")) {}
0020
0021 TriggerPrimitiveCollection DTBunchCrossingCleaner::clean(const TriggerPrimitiveCollection& inlist) const {
0022 TriggerPrimitiveCollection leftovers = inlist;
0023 TriggerPrimitiveCollection outlist;
0024
0025 auto tpin = inlist.cbegin();
0026 auto inend = inlist.cend();
0027 for (; tpin != inend; ++tpin) {
0028 const TriggerPrimitive::DTData data = tpin->getDTData();
0029
0030
0031 if (data.qualityCode != -1 && data.theta_quality != -1) {
0032 outlist.push_back(*tpin);
0033 auto toerase = std::find(leftovers.begin(), leftovers.end(), *tpin);
0034 if (toerase != leftovers.end()) {
0035 leftovers.erase(toerase);
0036 }
0037 }
0038
0039
0040
0041 if (data.qualityCode != -1 && data.theta_quality == -1) {
0042 auto tp_bx = leftovers.cbegin();
0043 auto tp_bx_end = leftovers.cend();
0044 for (; tp_bx != tp_bx_end; ++tp_bx) {
0045 if (*tp_bx == *tpin)
0046 continue;
0047 const TriggerPrimitive::DTData bx_data = tp_bx->getDTData();
0048
0049
0050 if (std::abs(bx_data.bx - data.bx) <= bx_window_size && bx_data.qualityCode == -1 &&
0051 bx_data.theta_quality != -1 && data.segment_number == bx_data.segment_number) {
0052
0053 L1MuDTChambPhDigi phi_digi(data.bx,
0054 data.wheel,
0055 data.sector,
0056 data.station,
0057 data.radialAngle,
0058 data.bendingAngle,
0059 data.qualityCode,
0060 data.Ts2TagCode,
0061 data.BxCntCode);
0062 int qual[7], position[7];
0063 for (int i = 0; i < 7; ++i) {
0064 qual[i] = 0;
0065 position[i] = 0;
0066 if (bx_data.theta_bti_group == i) {
0067 qual[i] = bx_data.theta_quality;
0068 position[i] = bx_data.segment_number;
0069 }
0070 }
0071 L1MuDTChambThDigi the_digi(data.bx, data.wheel, data.sector, data.station, position, qual);
0072
0073 DTChamberId the_id = tpin->detId<DTChamberId>();
0074 TriggerPrimitive newtp(the_id, phi_digi, the_digi, bx_data.theta_bti_group);
0075
0076 outlist.push_back(newtp);
0077
0078 auto phierase = std::find(leftovers.begin(), leftovers.end(), *tpin);
0079 auto theerase = std::find(leftovers.begin(), leftovers.end(), *tp_bx);
0080 if (phierase != leftovers.end()) {
0081 leftovers.erase(phierase);
0082 }
0083 if (theerase != leftovers.end()) {
0084 leftovers.erase(theerase);
0085 }
0086 break;
0087 }
0088 }
0089 }
0090 }
0091
0092 auto lo_tp = leftovers.cbegin();
0093 auto lo_end = leftovers.cend();
0094 for (; lo_tp != lo_end; ++lo_tp) {
0095 outlist.push_back(*lo_tp);
0096 }
0097 return outlist;
0098 }