CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes
HLTmmkFilter Class Reference

#include <HLTmmkFilter.h>

Inheritance diagram for HLTmmkFilter:
HLTFilter edm::global::EDFilter<> edm::global::EDFilterBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 HLTmmkFilter (const edm::ParameterSet &)
 
 ~HLTmmkFilter () override
 
- Public Member Functions inherited from HLTFilter
 HLTFilter (const edm::ParameterSet &config)
 
int module (edm::Event const &) const
 
const std::string * moduleLabel () const
 
int path (edm::Event const &) const
 
const std::string * pathName (edm::Event const &) const
 
std::pair< int, int > pmid (edm::Event const &) const
 
bool saveTags () const
 
 ~HLTFilter () override
 
- Public Member Functions inherited from edm::global::EDFilter<>
 EDFilter ()=default
 
bool hasAbilityToProduceInLumis () const final
 
bool hasAbilityToProduceInRuns () const final
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsStreamLuminosityBlocks () const final
 
bool wantsStreamRuns () const final
 
- Public Member Functions inherited from edm::global::EDFilterBase
 EDFilterBase ()
 
ModuleDescription const & moduleDescription () const
 
 ~EDFilterBase () override
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector< edm::ProductResolverIndex > const & indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
std::vector< edm::ProductResolverIndex > const & putTokenIndexToProductResolverIndex () const
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription const &)> registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, ModuleToResolverIndicies const &iIndicies, std::string const &moduleLabel)
 
 ~ProducerBase () noexcept(false) override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector< ProductResolverIndexAndSkipBit > const & itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 
- Static Public Member Functions inherited from HLTFilter
static void makeHLTFilterDescription (edm::ParameterSetDescription &desc)
 
- Static Public Member Functions inherited from edm::global::EDFilterBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 

Private Member Functions

void beginJob () override
 
void endJob () override
 
bool hltFilter (edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
 

Static Private Member Functions

static FreeTrajectoryState initialFreeState (const reco::Track &, const MagneticField *)
 
static int overlap (const reco::Candidate &, const reco::Candidate &)
 

Private Attributes

edm::InputTag beamSpotTag_
 
edm::EDGetTokenT< reco::BeamSpotbeamSpotToken_
 
const bool fastAccept_
 
const double maxEta_
 
const double maxInvMass_
 
const double maxNormalisedChi2_
 
const double minCosinePointingAngle_
 
const double minD0Significance_
 
const double minInvMass_
 
const double minLxySignificance_
 
const double minPt_
 
edm::InputTag muCandTag_
 
edm::EDGetTokenT< reco::RecoChargedCandidateCollectionmuCandToken_
 
const double thirdTrackMass_
 
edm::InputTag trkCandTag_
 
edm::EDGetTokenT< reco::RecoChargedCandidateCollectiontrkCandToken_
 

Additional Inherited Members

- Public Types inherited from edm::global::EDFilterBase
typedef EDFilterBase ModuleType
 
- Public Types inherited from edm::ProducerBase
using ModuleToResolverIndicies = std::unordered_multimap< std::string, std::tuple< edm::TypeID const *, const char *, edm::ProductResolverIndex >>
 
typedef ProductRegistryHelper::TypeLabelList TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

HLT Filter for b to (mumu) + X

Implementation: <Notes on="" implementation>="">

Definition at line 41 of file HLTmmkFilter.h.

Constructor & Destructor Documentation

HLTmmkFilter::HLTmmkFilter ( const edm::ParameterSet iConfig)
explicit

Definition at line 37 of file HLTmmkFilter.cc.

References ~HLTmmkFilter().

37  : HLTFilter(iConfig),
38  muCandTag_ (iConfig.getParameter<edm::InputTag>("MuCand")),
39  muCandToken_(consumes<reco::RecoChargedCandidateCollection>(muCandTag_)),
40  trkCandTag_ (iConfig.getParameter<edm::InputTag>("TrackCand")),
41  trkCandToken_(consumes<reco::RecoChargedCandidateCollection>(trkCandTag_)),
42  thirdTrackMass_(iConfig.getParameter<double>("ThirdTrackMass")),
43  maxEta_(iConfig.getParameter<double>("MaxEta")),
44  minPt_(iConfig.getParameter<double>("MinPt")),
45  minInvMass_(iConfig.getParameter<double>("MinInvMass")),
46  maxInvMass_(iConfig.getParameter<double>("MaxInvMass")),
47  maxNormalisedChi2_(iConfig.getParameter<double>("MaxNormalisedChi2")),
48  minLxySignificance_(iConfig.getParameter<double>("MinLxySignificance")),
49  minCosinePointingAngle_(iConfig.getParameter<double>("MinCosinePointingAngle")),
50  minD0Significance_(iConfig.getParameter<double>("MinD0Significance")),
51  fastAccept_(iConfig.getParameter<bool>("FastAccept")),
52  beamSpotTag_ (iConfig.getParameter<edm::InputTag> ("BeamSpotTag")),
53  beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_))
54 {
55  produces<VertexCollection>();
56  produces<CandidateCollection>();
57 }
T getParameter(std::string const &) const
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > trkCandToken_
Definition: HLTmmkFilter.h:59
const double maxInvMass_
Definition: HLTmmkFilter.h:65
edm::InputTag beamSpotTag_
Definition: HLTmmkFilter.h:71
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
Definition: HLTmmkFilter.h:72
edm::InputTag trkCandTag_
Definition: HLTmmkFilter.h:58
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > muCandToken_
Definition: HLTmmkFilter.h:57
HLTFilter(const edm::ParameterSet &config)
Definition: HLTFilter.cc:20
edm::InputTag muCandTag_
Definition: HLTmmkFilter.h:56
const double maxEta_
Definition: HLTmmkFilter.h:62
const double minCosinePointingAngle_
Definition: HLTmmkFilter.h:68
const double minInvMass_
Definition: HLTmmkFilter.h:64
const double thirdTrackMass_
Definition: HLTmmkFilter.h:61
const double minPt_
Definition: HLTmmkFilter.h:63
const double maxNormalisedChi2_
Definition: HLTmmkFilter.h:66
const double minD0Significance_
Definition: HLTmmkFilter.h:69
const bool fastAccept_
Definition: HLTmmkFilter.h:70
const double minLxySignificance_
Definition: HLTmmkFilter.h:67
HLTmmkFilter::~HLTmmkFilter ( )
overridedefault

Referenced by HLTmmkFilter().

Member Function Documentation

void HLTmmkFilter::beginJob ( void  )
overrideprivatevirtual

Reimplemented from edm::global::EDFilterBase.

Definition at line 84 of file HLTmmkFilter.cc.

84  {
85 
86 }
void HLTmmkFilter::endJob ( void  )
overrideprivatevirtual

Reimplemented from edm::global::EDFilterBase.

Definition at line 90 of file HLTmmkFilter.cc.

90  {
91 
92 }
void HLTmmkFilter::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 64 of file HLTmmkFilter.cc.

References edm::ConfigurationDescriptions::add(), edm::ParameterSetDescription::add(), and HLTFilter::makeHLTFilterDescription().

64  {
67  desc.add<edm::InputTag>("MuCand",edm::InputTag("hltMuTracks"));
68  desc.add<edm::InputTag>("TrackCand",edm::InputTag("hltMumukAllConeTracks"));
69  desc.add<double>("ThirdTrackMass",0.106);
70  desc.add<double>("MaxEta",2.5);
71  desc.add<double>("MinPt",3.0);
72  desc.add<double>("MinInvMass",1.2);
73  desc.add<double>("MaxInvMass",2.2);
74  desc.add<double>("MaxNormalisedChi2",10.0);
75  desc.add<double>("MinLxySignificance",3.0);
76  desc.add<double>("MinCosinePointingAngle",0.9);
77  desc.add<double>("MinD0Significance",0.0);
78  desc.add<bool>("FastAccept",false);
79  desc.add<edm::InputTag>("BeamSpotTag",edm::InputTag("hltOfflineBeamSpot"));
80  descriptions.add("hltmmkFilter",desc);
81 }
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void add(std::string const &label, ParameterSetDescription const &psetDescription)
bool HLTmmkFilter::hltFilter ( edm::Event iEvent,
const edm::EventSetup iSetup,
trigger::TriggerFilterObjectWithRefs filterproduct 
) const
overrideprivatevirtual

Implements HLTFilter.

Definition at line 96 of file HLTmmkFilter.cc.

References funct::abs(), accept(), trigger::TriggerFilterObjectWithRefs::addCollectionTag(), trigger::TriggerRefsCollections::addObject(), beamSpotToken_, counter, SoftLeptonByDistance_cfi::distance, reco::BeamSpot::dxdz(), reco::BeamSpot::dydz(), SiPixelPhase1TrackClustersV_cfi::e3, fastAccept_, edm::EventSetup::get(), edm::Ref< C, T, F >::get(), edm::Event::getByToken(), trigger::TriggerRefsCollections::getObjects(), mps_fire::i, initialFreeState(), TransientVertex::isValid(), LogDebug, maxEta_, maxInvMass_, maxNormalisedChi2_, minCosinePointingAngle_, minD0Significance_, minInvMass_, minLxySignificance_, minPt_, eostools::move(), muCandTag_, muCandToken_, TransientVertex::normalisedChiSquared(), convertSQLitetoXML_cfg::output, overlap(), AlCaHLTBitMon_ParallelJobs::p, p1, p2, p3, PV3DBase< T, PVType, FrameType >::perp(), TransientVertex::position(), TransientVertex::positionError(), edm::ESHandle< T >::product(), edm::Event::put(), GlobalErrorBase< T, ErrorWeightType >::rerr(), HLTFilter::saveTags(), Measurement1D::significance(), createPayload::skip, mathSSE::sqrt(), thirdTrackMass_, TrajectoryStateClosestToBeamLine::transverseImpactParameter(), trigger::TriggerMuon, trigger::TriggerTrack, trkCandTag_, trkCandToken_, KalmanVertexFitter::vertex(), particleFlowSuperClusterECAL_cfi::vertexCollection, PV3DBase< T, PVType, FrameType >::x(), reco::BeamSpot::x0(), PV3DBase< T, PVType, FrameType >::y(), reco::BeamSpot::y0(), PV3DBase< T, PVType, FrameType >::z(), and reco::BeamSpot::z0().

96  {
97 
98  const double MuMass(0.106);
99  const double MuMass2(MuMass*MuMass);
100 
101  const double thirdTrackMass2(thirdTrackMass_*thirdTrackMass_);
102 
103  unique_ptr<CandidateCollection> output(new CandidateCollection());
104  unique_ptr<VertexCollection> vertexCollection(new VertexCollection());
105 
106  //get the transient track builder:
108  iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder",theB);
109 
110  //get the beamspot position
111  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
112  iEvent.getByToken(beamSpotToken_,recoBeamSpotHandle);
113  const reco::BeamSpot& vertexBeamSpot = *recoBeamSpotHandle;
114 
115  ESHandle<MagneticField> bFieldHandle;
116  iSetup.get<IdealMagneticFieldRecord>().get(bFieldHandle);
117 
118  const MagneticField* magField = bFieldHandle.product();
119 
120  TSCBLBuilderNoMaterial blsBuilder;
121 
122  // Ref to Candidate object to be recorded in filter object
125  RecoChargedCandidateRef refTrk;
126 
127  // get hold of muon trks
129  iEvent.getByToken(muCandToken_,mucands);
130 
131  // get track candidates around displaced muons
133  iEvent.getByToken(trkCandToken_,trkcands);
134 
135  if (saveTags()) {
136  filterproduct.addCollectionTag(muCandTag_);
137  filterproduct.addCollectionTag(trkCandTag_);
138  }
139 
140  double e1,e2,e3;
142 
143  //TrackRefs to mu cands in trkcand
144  vector<TrackRef> trkMuCands;
145 
146  //Already used mu tracks to avoid double counting candidates
147  vector<bool> isUsedCand(trkcands->size(),false);
148 
149  int counter = 0;
150 
151  //run algorithm
152  for (auto mucand1=mucands->begin(), endCand1=mucands->end(); mucand1!=endCand1; ++mucand1) {
153 
154  if ( mucands->size()<2) break;
155  if ( trkcands->empty()) break;
156 
157  TrackRef trk1 = mucand1->get<TrackRef>();
158  LogDebug("HLTDisplacedMumukFilter") << " 1st muon: q*pt= " << trk1->charge()*trk1->pt() << ", eta= " << trk1->eta() << ", hits= " << trk1->numberOfValidHits();
159 
160  // eta cut
161  if (fabs(trk1->eta()) > maxEta_) continue;
162 
163  // Pt threshold cut
164  if (trk1->pt() < minPt_) continue;
165 
166  auto mucand2 = mucand1; ++mucand2;
167 
168  for (auto endCand2=mucands->end(); mucand2!=endCand2; ++mucand2) {
169 
170  TrackRef trk2 = mucand2->get<TrackRef>();
171 
172  LogDebug("HLTDisplacedMumukFilter") << " 2nd muon: q*pt= " << trk2->charge()*trk2->pt() << ", eta= " << trk2->eta() << ", hits= " << trk2->numberOfValidHits();
173 
174  // eta cut
175  if (fabs(trk2->eta()) > maxEta_) continue;
176 
177  // Pt threshold cut
178  if (trk2->pt() < minPt_) continue;
179 
180  RecoChargedCandidateCollection::const_iterator trkcand, endCandTrk;
181 
182  std::vector<bool>::iterator isUsedIter, endIsUsedCand;
183 
184  //get overlapping muon candidates
185  for ( trkcand = trkcands->begin(), endCandTrk=trkcands->end(), isUsedIter = isUsedCand.begin(), endIsUsedCand = isUsedCand.end(); trkcand != endCandTrk && isUsedIter != endIsUsedCand; ++trkcand, ++isUsedIter) {
186  TrackRef trk3 = trkcand->get<TrackRef>();
187 
188  //check for overlapping muon tracks and save TrackRefs
189  if (overlap(*mucand1,*trkcand)) {
190  trkMuCands.push_back(trk3);
191  *isUsedIter = true;
192  continue;
193  }
194  else if (overlap(*mucand2,*trkcand)){
195  trkMuCands.push_back(trk3);
196  *isUsedIter = true;
197  continue;
198  }
199 
200  if(trkMuCands.size()==2) break;
201  }
202 
203  //Not all overlapping candidates found, skip event
204  //if (trkMuCands.size()!=2) continue;
205 
206  //combine muons with all tracks
207  for ( trkcand = trkcands->begin(), endCandTrk=trkcands->end(), isUsedIter = isUsedCand.begin(), endIsUsedCand = isUsedCand.end(); trkcand != endCandTrk && isUsedIter != endIsUsedCand; ++trkcand, ++isUsedIter) {
208 
209  TrackRef trk3 = trkcand->get<TrackRef>();
210 
211  LogDebug("HLTDisplacedMumukFilter") << " 3rd track: q*pt= " << trk3->charge()*trk3->pt() << ", eta= " << trk3->eta() << ", hits= " << trk3->numberOfValidHits();
212 
213  //skip overlapping muon candidates
214  bool skip=false;
215  for (auto & trkMuCand : trkMuCands) if(trk3==trkMuCand) skip=true;
216  if(skip) continue;
217 
218  //skip already used tracks
219  if(*isUsedIter) continue;
220 
221  // eta cut
222  if (fabs(trk3->eta()) > maxEta_) continue;
223 
224  // Pt threshold cut
225  if (trk3->pt() < minPt_) continue;
226 
227  // Combined system
228  e1 = sqrt(trk1->momentum().Mag2()+MuMass2);
229  e2 = sqrt(trk2->momentum().Mag2()+MuMass2);
230  e3 = sqrt(trk3->momentum().Mag2()+thirdTrackMass2);
231 
232  p1 = Particle::LorentzVector(trk1->px(),trk1->py(),trk1->pz(),e1);
233  p2 = Particle::LorentzVector(trk2->px(),trk2->py(),trk2->pz(),e2);
234  p3 = Particle::LorentzVector(trk3->px(),trk3->py(),trk3->pz(),e3);
235 
236  p = p1+p2+p3;
237 
238  //invariant mass cut
239  double invmass = abs(p.mass());
240 
241  LogDebug("HLTDisplacedMumukFilter") << " Invmass= " << invmass;
242 
243  if (invmass>maxInvMass_ || invmass<minInvMass_) continue;
244 
245  // do the vertex fit
246  vector<TransientTrack> t_tks;
247  t_tks.push_back((*theB).build(&trk1));
248  t_tks.push_back((*theB).build(&trk2));
249  t_tks.push_back((*theB).build(&trk3));
250 
251  if (t_tks.size()!=3) continue;
252 
253  FreeTrajectoryState InitialFTS = initialFreeState(*trk3, magField);
254  TrajectoryStateClosestToBeamLine tscb( blsBuilder(InitialFTS, *recoBeamSpotHandle) );
255  double d0sig = tscb.transverseImpactParameter().significance();
256 
257  if (d0sig < minD0Significance_) continue;
258 
259  KalmanVertexFitter kvf;
260  TransientVertex tv = kvf.vertex(t_tks);
261 
262  if (!tv.isValid()) continue;
263 
264  Vertex vertex = tv;
265 
266  // get vertex position and error to calculate the decay length significance
267  GlobalPoint secondaryVertex = tv.position();
268  GlobalError err = tv.positionError();
269 
270  //calculate decay length significance w.r.t. the beamspot
271  GlobalPoint displacementFromBeamspot( -1*((vertexBeamSpot.x0() -secondaryVertex.x()) + (secondaryVertex.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dxdz()), -1*((vertexBeamSpot.y0() - secondaryVertex.y())+ (secondaryVertex.z() -vertexBeamSpot.z0()) * vertexBeamSpot.dydz()), 0);
272 
273  float lxy = displacementFromBeamspot.perp();
274  float lxyerr = sqrt(err.rerr(displacementFromBeamspot));
275 
276  // get normalizes chi2
277  float normChi2 = tv.normalisedChiSquared();
278 
279  //calculate the angle between the decay length and the mumu momentum
280  Vertex::Point vperp(displacementFromBeamspot.x(),displacementFromBeamspot.y(),0.);
281  math::XYZVector pperp(p.x(),p.y(),0.);
282 
283  float cosAlpha = vperp.Dot(pperp)/(vperp.R()*pperp.R());
284 
285  LogDebug("HLTDisplacedMumukFilter") << " vertex fit normalised chi2: " << normChi2 << ", Lxy significance: " << lxy/lxyerr << ", cosine pointing angle: " << cosAlpha;
286 
287  // put vertex in the event
288  vertexCollection->push_back(vertex);
289 
290  if (normChi2 > maxNormalisedChi2_) continue;
291  if (lxy/lxyerr < minLxySignificance_) continue;
292  if(cosAlpha < minCosinePointingAngle_) continue;
293 
294  LogDebug("HLTDisplacedMumukFilter") << " Event passed!";
295 
296  //Add event
297  ++counter;
298 
299  //Get refs
300  bool i1done = false;
301  bool i2done = false;
302  bool i3done = false;
303  vector<RecoChargedCandidateRef> vref;
304  filterproduct.getObjects(TriggerMuon,vref);
305  for (auto & i : vref) {
307  TrackRef trktmp = candref->get<TrackRef>();
308  if (trktmp==trk1) {
309  i1done = true;
310  } else if (trktmp==trk2) {
311  i2done = true;
312  } else if (trktmp==trk3) {
313  i3done = true;
314  }
315  if (i1done && i2done && i3done) break;
316  }
317 
318  if (!i1done) {
319  refMu1=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (mucands,distance(mucands->begin(), mucand1)));
320  filterproduct.addObject(TriggerMuon,refMu1);
321  }
322  if (!i2done) {
323  refMu2=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (mucands,distance(mucands->begin(),mucand2)));
324  filterproduct.addObject(TriggerMuon,refMu2);
325  }
326  if (!i3done) {
327  refTrk=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (trkcands,distance(trkcands->begin(),trkcand)));
328  filterproduct.addObject(TriggerTrack,refTrk);
329  }
330 
331  if (fastAccept_) break;
332  }
333 
334  trkMuCands.clear();
335  }
336  }
337 
338  // filter decision
339  const bool accept (counter >= 1);
340 
341  LogDebug("HLTDisplacedMumukFilter") << " >>>>> Result of HLTDisplacedMumukFilter is "<< accept << ", number of muon pairs passing thresholds= " << counter;
342 
343  iEvent.put(std::move(vertexCollection));
344 
345  return accept;
346 }
#define LogDebug(id)
GlobalError positionError() const
static int overlap(const reco::Candidate &, const reco::Candidate &)
double z0() const
z coordinate
Definition: BeamSpot.h:68
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > trkCandToken_
Definition: HLTmmkFilter.h:59
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:125
void getObjects(Vids &ids, VRphoton &refs) const
various physics-level getters:
const double maxInvMass_
Definition: HLTmmkFilter.h:65
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
T y() const
Definition: PV3DBase.h:63
edm::Ref< RecoChargedCandidateCollection > RecoChargedCandidateRef
reference to an object in a collection of RecoChargedCandidate objects
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
double dydz() const
dydz slope
Definition: BeamSpot.h:84
float normalisedChiSquared() const
T sqrt(T t)
Definition: SSEVec.h:18
GlobalPoint position() const
T z() const
Definition: PV3DBase.h:64
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
CachingVertex< 5 > vertex(const std::vector< reco::TransientTrack > &tracks) const override
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
Definition: HLTmmkFilter.h:72
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:243
edm::InputTag trkCandTag_
Definition: HLTmmkFilter.h:58
double p2[4]
Definition: TauolaWrapper.h:90
double dxdz() const
dxdz slope
Definition: BeamSpot.h:82
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > muCandToken_
Definition: HLTmmkFilter.h:57
edm::InputTag muCandTag_
Definition: HLTmmkFilter.h:56
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
const double maxEta_
Definition: HLTmmkFilter.h:62
T rerr(const GlobalPoint &aPoint) const
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
const double minCosinePointingAngle_
Definition: HLTmmkFilter.h:68
const double minInvMass_
Definition: HLTmmkFilter.h:64
const double thirdTrackMass_
Definition: HLTmmkFilter.h:61
bool saveTags() const
Definition: HLTFilter.h:45
double p1[4]
Definition: TauolaWrapper.h:89
static std::atomic< unsigned int > counter
T get() const
Definition: EventSetup.h:71
double y0() const
y coordinate
Definition: BeamSpot.h:66
const double minPt_
Definition: HLTmmkFilter.h:63
const double maxNormalisedChi2_
Definition: HLTmmkFilter.h:66
static FreeTrajectoryState initialFreeState(const reco::Track &, const MagneticField *)
edm::OwnVector< Candidate > CandidateCollection
collection of Candidate objects
Definition: CandidateFwd.h:21
T x() const
Definition: PV3DBase.h:62
const double minD0Significance_
Definition: HLTmmkFilter.h:69
T const * product() const
Definition: ESHandle.h:86
bool isValid() const
def move(src, dest)
Definition: eostools.py:511
math::PtEtaPhiELorentzVectorF LorentzVector
const bool fastAccept_
Definition: HLTmmkFilter.h:70
const double minLxySignificance_
Definition: HLTmmkFilter.h:67
double p3[4]
Definition: TauolaWrapper.h:91
double x0() const
x coordinate
Definition: BeamSpot.h:64
FreeTrajectoryState HLTmmkFilter::initialFreeState ( const reco::Track tk,
const MagneticField field 
)
staticprivate

Definition at line 349 of file HLTmmkFilter.cc.

References reco::TrackBase::charge(), reco::TrackBase::covariance(), reco::TrackBase::momentum(), and reco::TrackBase::vertex().

Referenced by hltFilter().

351 {
353  GlobalPoint gpos( pos);
354  Basic3DVector<float> mom( tk.momentum());
355  GlobalVector gmom( mom);
356  GlobalTrajectoryParameters par( gpos, gmom, tk.charge(), field);
358  return FreeTrajectoryState( par, err);
359 }
const Vector & momentum() const
track momentum vector
Definition: TrackBase.h:714
const Point & vertex() const
reference point on the track. This method is DEPRECATED, please use referencePoint() instead ...
Definition: TrackBase.h:738
CovarianceMatrix covariance() const
return track covariance matrix
Definition: TrackBase.h:782
int charge() const
track electric charge
Definition: TrackBase.h:606
int HLTmmkFilter::overlap ( const reco::Candidate a,
const reco::Candidate b 
)
staticprivate

Definition at line 361 of file HLTmmkFilter.cc.

References hiPixelPairStep_cff::deltaPhi, MillePedeFileConverter_cfg::e, reco::Candidate::eta(), reco::Candidate::phi(), and reco::Candidate::pt().

Referenced by hltFilter().

361  {
362 
363  double eps(1.44e-4);
364 
365  double dpt = a.pt() - b.pt();
366  dpt *= dpt;
367 
368  double dphi = deltaPhi(a.phi(), b.phi());
369  dphi *= dphi;
370 
371  double deta = a.eta() - b.eta();
372  deta *= deta;
373 
374  if ((dphi + deta) < eps) {
375  return 1;
376  }
377 
378  return 0;
379 
380 }
virtual double eta() const =0
momentum pseudorapidity
virtual double pt() const =0
transverse momentum
virtual double phi() const =0
momentum azimuthal angle

Member Data Documentation

edm::InputTag HLTmmkFilter::beamSpotTag_
private

Definition at line 71 of file HLTmmkFilter.h.

edm::EDGetTokenT<reco::BeamSpot> HLTmmkFilter::beamSpotToken_
private

Definition at line 72 of file HLTmmkFilter.h.

Referenced by hltFilter().

const bool HLTmmkFilter::fastAccept_
private

Definition at line 70 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::maxEta_
private

Definition at line 62 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::maxInvMass_
private

Definition at line 65 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::maxNormalisedChi2_
private

Definition at line 66 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::minCosinePointingAngle_
private

Definition at line 68 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::minD0Significance_
private

Definition at line 69 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::minInvMass_
private

Definition at line 64 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::minLxySignificance_
private

Definition at line 67 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::minPt_
private

Definition at line 63 of file HLTmmkFilter.h.

Referenced by hltFilter().

edm::InputTag HLTmmkFilter::muCandTag_
private

Definition at line 56 of file HLTmmkFilter.h.

Referenced by hltFilter().

edm::EDGetTokenT<reco::RecoChargedCandidateCollection> HLTmmkFilter::muCandToken_
private

Definition at line 57 of file HLTmmkFilter.h.

Referenced by hltFilter().

const double HLTmmkFilter::thirdTrackMass_
private

Definition at line 61 of file HLTmmkFilter.h.

Referenced by hltFilter().

edm::InputTag HLTmmkFilter::trkCandTag_
private

Definition at line 58 of file HLTmmkFilter.h.

Referenced by hltFilter().

edm::EDGetTokenT<reco::RecoChargedCandidateCollection> HLTmmkFilter::trkCandToken_
private

Definition at line 59 of file HLTmmkFilter.h.

Referenced by hltFilter().