CMS 3D CMS Logo

HLTCSCRing2or3Filter.cc
Go to the documentation of this file.
1 #include "HLTCSCRing2or3Filter.h"
2 
6 
12 
15 
17  , m_input(iConfig.getParameter<edm::InputTag>("input"))
18  , m_minHits(iConfig.getParameter<unsigned int>("minHits"))
19  , m_xWindow(iConfig.getParameter<double>("xWindow"))
20  , m_yWindow(iConfig.getParameter<double>("yWindow"))
21 {
22  cscrechitsToken = consumes<CSCRecHit2DCollection>(m_input);
23 }
24 
26 
27 void
31  desc.add<edm::InputTag>("input",edm::InputTag("hltCsc2DRecHits"));
32  desc.add<unsigned int>("minHits",4);
33  desc.add<double>("xWindow",2.);
34  desc.add<double>("yWindow",2.);
35  descriptions.add("hltCSCRing2or3Filter",desc);
36 }
37 
40  iEvent.getByToken(cscrechitsToken, hits);
41 
42  edm::ESHandle<CSCGeometry> cscGeometry;
43  bool got_cscGeometry = false;
44 
45  std::map<int, std::vector<const CSCRecHit2D*> > chamber_tohit;
46 
47  for (auto const & hit : *hits) {
48  CSCDetId id(hit.geographicalId());
49  int chamber_id = CSCDetId(id.endcap(), id.station(), id.ring(), id.chamber(), 0).rawId();
50 
51  if (id.ring() == 2 || id.ring() == 3) {
52  std::map<int, std::vector<const CSCRecHit2D*> >::const_iterator chamber_iter = chamber_tohit.find(chamber_id);
53  if (chamber_iter == chamber_tohit.end()) {
54  std::vector<const CSCRecHit2D*> newlist;
55  newlist.push_back(&hit);
56  }
57  chamber_tohit[chamber_id].push_back(&hit);
58  } // end if this ring is selected
59  } // end loop over hits
60 
61  unsigned int minHitsAlmostSquared = (m_minHits-1) * (m_minHits-2);
62  for (std::map<int, std::vector<const CSCRecHit2D*> >::const_iterator chamber_iter = chamber_tohit.begin();
63  chamber_iter != chamber_tohit.end();
64  ++chamber_iter) {
65 
66  if (chamber_iter->second.size() >= m_minHits) {
67  if (!got_cscGeometry) {
68  iSetup.get<MuonGeometryRecord>().get(cscGeometry);
69  got_cscGeometry = true;
70  }
71 
72  unsigned int pairs_in_window = 0;
73  for (auto hit1 = chamber_iter->second.begin(); hit1 != chamber_iter->second.end(); ++hit1) {
74  for (auto hit2 = chamber_iter->second.begin(); hit2 != hit1; ++hit2) {
75  GlobalPoint pos1 = cscGeometry->idToDet((*hit1)->geographicalId())->surface().toGlobal((*hit1)->localPosition());
76  GlobalPoint pos2 = cscGeometry->idToDet((*hit2)->geographicalId())->surface().toGlobal((*hit2)->localPosition());
77 
78  if (fabs(pos1.x() - pos2.x()) < m_xWindow && fabs(pos1.y() - pos2.y()) < m_yWindow) pairs_in_window++;
79 
80  if (pairs_in_window >= minHitsAlmostSquared) return true;
81  } // end loop over hits
82  } // end loop over hits
83 
84  } // end if chamber has enough hits
85  } // end loop over chambers
86 
87  return false;
88 }
89 
90 // declare this class as a framework plugin
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
const GeomDet * idToDet(DetId) const override
Definition: CSCGeometry.cc:114
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:54
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:50
T y() const
Definition: PV3DBase.h:63
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void add(std::string const &label, ParameterSetDescription const &psetDescription)
~HLTCSCRing2or3Filter() override
HLT enums.
edm::EDGetTokenT< CSCRecHit2DCollection > cscrechitsToken
T get() const
Definition: EventSetup.h:71
HLTCSCRing2or3Filter(const edm::ParameterSet &)
T x() const
Definition: PV3DBase.h:62