CMS 3D CMS Logo

HLTMuonL2FromL1TPreFilter.cc
Go to the documentation of this file.
1 
21 
22 //
23 // constructors and destructor
24 //
25 
27  beamSpotTag_( iConfig.getParameter<edm::InputTag>("BeamSpotTag") ),
28  beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_)),
29  candTag_( iConfig.getParameter<edm::InputTag >("CandTag") ),
30  candToken_(consumes<reco::RecoChargedCandidateCollection>(candTag_)),
31  previousCandTag_( iConfig.getParameter<edm::InputTag >("PreviousCandTag") ),
32  previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
33  seedMapTag_( iConfig.getParameter<edm::InputTag >("SeedMapTag") ),
34  seedMapToken_(consumes<SeedMap>(seedMapTag_)),
35  minN_( iConfig.getParameter<int>("MinN") ),
36  maxEta_( iConfig.getParameter<double>("MaxEta") ),
37  absetaBins_( iConfig.getParameter<std::vector<double> >("AbsEtaBins") ),
38  minNstations_( iConfig.getParameter<std::vector<int> >("MinNstations") ),
39  minNhits_( iConfig.getParameter<std::vector<int> >("MinNhits") ),
40  cutOnChambers_( iConfig.getParameter<bool>("CutOnChambers") ),
41  minNchambers_( iConfig.getParameter<std::vector<int> >("MinNchambers") ),
42  maxDr_( iConfig.getParameter<double>("MaxDr") ),
43  minDr_( iConfig.getParameter<double>("MinDr") ),
44  maxDz_( iConfig.getParameter<double>("MaxDz") ),
45  min_DxySig_(iConfig.getParameter<double> ("MinDxySig")),
46  minPt_( iConfig.getParameter<double>("MinPt") ),
47  nSigmaPt_( iConfig.getParameter<double>("NSigmaPt") ),
48  matchPreviousCand_( iConfig.getParameter<bool>("MatchToPreviousCand") )
49 {
50  using namespace std;
51 
52  // check that number of eta bins matches number of nStation cuts
53  if( minNstations_.size()!=absetaBins_.size() ||
54  minNhits_.size()!=absetaBins_.size() ||
55  ( cutOnChambers_ && minNchambers_.size()!=absetaBins_.size() ) ) {
56  throw cms::Exception("Configuration") << "Number of MinNstations, MinNhits, or MinNchambers cuts "
57  << "does not match number of eta bins." << endl;
58  }
59 
60  if(absetaBins_.size()>1) {
61  for(unsigned int i=0; i<absetaBins_.size()-1; ++i) {
62  if(absetaBins_[i+1]<=absetaBins_[i])
63  throw cms::Exception("Configuration") << "Absolute eta bins must be in increasing order." << endl;
64  }
65  }
66 
67  // dump parameters for debugging
68  if(edm::isDebugEnabled()){
69  ostringstream ss;
70  ss<<"Constructed with parameters:"<<endl;
71  ss<<" BeamSpotTag = "<<beamSpotTag_.encode()<<endl;
72  ss<<" CandTag = "<<candTag_.encode()<<endl;
73  ss<<" PreviousCandTag = "<<previousCandTag_.encode()<<endl;
74  ss<<" SeedMapTag = "<<seedMapTag_.encode()<<endl;
75  ss<<" MinN = "<<minN_<<endl;
76  ss<<" MaxEta = "<<maxEta_<<endl;
77  ss<<" MinNstations = ";
78  for(unsigned int j=0; j<absetaBins_.size(); ++j) {
79  ss<<minNstations_[j]<<" (|eta|<"<<absetaBins_[j]<<"), ";
80  }
81  ss<<endl;
82  ss<<" MinNhits = ";
83  for(unsigned int j=0; j<absetaBins_.size(); ++j) {
84  ss<<minNhits_[j]<<" (|eta|<"<<absetaBins_[j]<<"), ";
85  }
86  ss<<endl;
87  ss<<" CutOnChambers = " <<cutOnChambers_<<endl;
88  if ( cutOnChambers_ ) {
89  ss<<" MinNchambers = ";
90  for(unsigned int j=0; j<absetaBins_.size(); ++j) {
91  ss<<minNchambers_[j]<<" (|eta|<"<<absetaBins_[j]<<"), ";
92  }
93  }
94  ss<<endl;
95  ss<<" MaxDr = "<<maxDr_<<endl;
96  ss<<" MinDr = "<<minDr_<<endl;
97  ss<<" MaxDz = "<<maxDz_<<endl;
98  ss<<" MinDxySig = "<<min_DxySig_<<endl;
99  ss<<" MinPt = "<<minPt_<<endl;
100  ss<<" NSigmaPt = "<<nSigmaPt_<<endl;
101  ss<<" saveTags= "<<saveTags();
102  LogDebug("HLTMuonL2FromL1TPreFilter")<<ss.str();
103  }
104 }
105 
107 
108 void
112  desc.add<edm::InputTag>("BeamSpotTag",edm::InputTag("hltOfflineBeamSpot"));
113  desc.add<edm::InputTag>("CandTag",edm::InputTag("hltL2MuonCandidates"));
114  desc.add<edm::InputTag>("PreviousCandTag",edm::InputTag(""));
115  desc.add<edm::InputTag>("SeedMapTag",edm::InputTag("hltL2Muons"));
116  desc.add<int>("MinN",1);
117  desc.add<double>("MaxEta",2.5);
118  desc.add<std::vector<double> >("AbsEtaBins", std::vector<double>(1, 9999.));
119  desc.add<std::vector<int> >("MinNstations", std::vector<int>(1, 1));
120  desc.add<std::vector<int> >("MinNhits", std::vector<int>(1, 0));
121  desc.add<bool> ("CutOnChambers", false);
122  desc.add<std::vector<int> >("MinNchambers", std::vector<int>(1, 0));
123  desc.add<double>("MaxDr",9999.0);
124  desc.add<double>("MinDr",-1.0);
125  desc.add<double>("MaxDz",9999.0);
126  desc.add<double>("MinDxySig",-1.0);
127  desc.add<double>("MinPt",0.0);
128  desc.add<double>("NSigmaPt",0.0);
129  desc.add<bool>("MatchToPreviousCand",true);
130  descriptions.add("hltMuonL2FromL1TPreFilter",desc);
131 }
132 
133 //
134 // member functions
135 //
136 
137 // ------------ method called to produce the data ------------
139 {
140  // All HLT filters must create and fill an HLT filter object,
141  // recording any reconstructed physics objects satisfying (or not)
142  // this HLT filter, and place it in the Event.
143 
144  using namespace std;
145  using namespace edm;
146  using namespace reco;
147  using namespace trigger;
148 // using namespace l1extra;
149 
150  // save Tag
151  if (saveTags()) filterproduct.addCollectionTag(candTag_);
152 
153  // get hold of all muon candidates available at this level
155  iEvent.getByToken(candToken_, allMuons);
156 
157  // get hold of the beam spot
158  Handle<BeamSpot> beamSpotHandle;
159  iEvent.getByToken(beamSpotToken_, beamSpotHandle);
160  BeamSpot::Point beamSpot = beamSpotHandle->position();
161 
162  // get the L2 to L1 map object for this event
164 
165  // number of eta bins for cut on number of stations
166  const std::vector<double>::size_type nAbsetaBins = absetaBins_.size();
167 
168  // look at all allMuons, check cuts and add to filter object
169  int n = 0;
170  for(auto cand=allMuons->begin(); cand!=allMuons->end(); cand++){
171  TrackRef mu = cand->get<TrackRef>();
172 
173  // check if this muon passed previous level
174  if(matchPreviousCand_ && !mapL2ToL1.isTriggeredByL1(mu)) continue;
175 
176  // eta cut
177  if(std::abs(mu->eta()) > maxEta_) continue;
178 
179  // cut on number of stations
180  bool failNstations(false), failNhits(false), failNchambers(false);
181  for(unsigned int i=0; i<nAbsetaBins; ++i) {
182  if( std::abs(mu->eta())<absetaBins_[i] ) {
183  if(mu->hitPattern().muonStationsWithAnyHits() < minNstations_[i]) {
184  failNstations=true;
185  }
186  if(mu->numberOfValidHits() < minNhits_[i]) {
187  failNhits=true;
188  }
189  if( cutOnChambers_ &&
190  ( mu->hitPattern().dtStationsWithAnyHits() +
191  mu->hitPattern().cscStationsWithAnyHits() < minNchambers_[i]) ) {
192  failNchambers=true;
193  }
194  break;
195  }
196  }
197  if(failNstations || failNhits || failNchambers) continue;
198 
199  //dr cut
200  if(std::abs(mu->dxy(beamSpot)) > maxDr_) continue;
201 
202  //dr cut
203  if(std::abs(mu->dxy(beamSpot)) < minDr_) continue;
204 
205  //dz cut
206  if(std::abs(mu->dz(beamSpot)) > maxDz_) continue;
207 
208  // dxy significance cut (safeguard against bizarre values)
209  if (min_DxySig_ > 0 && (mu->dxyError() <= 0 || std::abs(mu->dxy(beamSpot)/mu->dxyError()) < min_DxySig_)) continue;
210 
211  // Pt threshold cut
212  double pt = mu->pt();
213  double abspar0 = std::abs(mu->parameter(0));
214  double ptLx = pt;
215  // convert 50% efficiency threshold to 90% efficiency threshold
216  if(abspar0 > 0) ptLx += nSigmaPt_*mu->error(0)/abspar0*pt;
217  if(ptLx < minPt_) continue;
218 
219  // add the good candidate to the filter object
220  filterproduct.addObject(TriggerMuon, RecoChargedCandidateRef(Ref<RecoChargedCandidateCollection>(allMuons, cand-allMuons->begin())));
221 
222  n++;
223  }
224 
225  // filter decision
226  const bool accept (n >= minN_);
227 
228  // dump event for debugging
229  if(edm::isDebugEnabled()){
230  ostringstream ss;
231  ss<<"L2mu#"
232  <<'\t'<<"q*pt"<<'\t' //scientific is too wide
233  <<'\t'<<"q*ptLx"<<'\t' //scientific is too wide
234  <<'\t'<<"eta"
235  <<'\t'<<"phi"
236  <<'\t'<<"nStations"
237  <<'\t'<<"nHits"
238  <<'\t'<<"dr"<<'\t' //scientific is too wide
239  <<'\t'<<"dz"<<'\t' //scientific is too wide
240  <<'\t'<<"L1seed#"
241  <<'\t'<<"isPrev"
242  <<'\t'<<"isFired"
243  <<endl;
244  ss<<"-----------------------------------------------------------------------------------------------------------------------"<<endl;
245  for (auto cand = allMuons->begin(); cand != allMuons->end(); cand++) {
246  TrackRef mu = cand->get<TrackRef>();
247  ss<<setprecision(2)
248  <<cand-allMuons->begin()
249  <<'\t'<<scientific<<mu->charge()*mu->pt()
250  <<'\t'<<scientific<<mu->charge()*mu->pt()*(1. + ((mu->parameter(0) != 0) ? nSigmaPt_*mu->error(0)/std::abs(mu->parameter(0)) : 0.))
251  <<'\t'<<fixed<<mu->eta()
252  <<'\t'<<fixed<<mu->phi()
253  <<'\t'<<mu->hitPattern().muonStationsWithAnyHits()
254  <<'\t'<<mu->numberOfValidHits()
255  <<'\t'<<scientific<<mu->d0()
256  <<'\t'<<scientific<<mu->dz()
257  <<'\t'<<mapL2ToL1.getL1Keys(mu)
258  <<'\t'<<mapL2ToL1.isTriggeredByL1(mu);
259  vector<RecoChargedCandidateRef> firedMuons;
260  filterproduct.getObjects(TriggerMuon, firedMuons);
261  ss<<'\t'<<(find(firedMuons.begin(), firedMuons.end(), RecoChargedCandidateRef(Ref<RecoChargedCandidateCollection>(allMuons, cand-allMuons->begin()))) != firedMuons.end())
262  <<endl;
263  }
264  ss<<"-----------------------------------------------------------------------------------------------------------------------"<<endl;
265  ss<<"Decision of filter is "<<accept<<", number of muons passing = "<<filterproduct.muonSize();
266  LogDebug("HLTMuonL2FromL1TPreFilter")<<ss.str();
267  }
268 
269  return accept;
270 }
271 
272 
273 // declare this class as a framework plugin
#define LogDebug(id)
bool isDebugEnabled()
edm::InputTag previousCandTag_
input tag of the preceeding L1 filter in the path
void getObjects(Vids &ids, VRphoton &refs) const
various physics-level getters:
~HLTMuonL2FromL1TPreFilter() override
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > candToken_
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
bool isTriggeredByL1(reco::TrackRef &l2muon)
checks if a L2 muon was seeded by a fired L1
std::vector< int > minNstations_
minimum number of muon stations used
double minDr_
cut on impact parameter wrt to the beam spot
edm::Ref< RecoChargedCandidateCollection > RecoChargedCandidateRef
reference to an object in a collection of RecoChargedCandidate objects
math::XYZPoint Point
point in the space
Definition: BeamSpot.h:29
bool cutOnChambers_
choose whether to apply cut on number of chambers (DT+CSC)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
uint16_t size_type
std::string encode() const
Definition: InputTag.cc:159
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
bool matchPreviousCand_
require the matching with the L1 firing the L1 filter
edm::EDGetTokenT< SeedMap > seedMapToken_
double min_DxySig_
dxy significance cut
std::vector< int > minNhits_
minimum number of valid muon hits
edm::InputTag beamSpotTag_
input tag of the beam spot
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
const int mu
Definition: Constants.h:22
ParameterDescriptionBase * add(U const &iLabel, T const &value)
int minN_
minimum number of muons to fire the trigger
double minPt_
pt threshold in GeV
std::vector< RecoChargedCandidate > RecoChargedCandidateCollection
collectin of RecoChargedCandidate objects
edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > previousCandToken_
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
std::string getL1Keys(reco::TrackRef &l2muon)
returns the indices of L1 seeds
std::vector< int > minNchambers_
minimum number of valid chambers
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
double maxDz_
cut on dz wrt to the beam spot
HLTMuonL2FromL1TPreFilter(const edm::ParameterSet &)
double maxDr_
cut on impact parameter wrt to the beam spot
const Point & position() const
position
Definition: BeamSpot.h:62
edm::InputTag candTag_
input tag of L2 muons
double nSigmaPt_
pt uncertainty margin (in number of sigmas)
edm::InputTag seedMapTag_
input tag of the map from the L2 seed to the sister L2 seeds of cleaned tracks