CMS 3D CMS Logo

HLTDisplacedmumuFilter.cc
Go to the documentation of this file.
1 #include <iostream>
2 
5 
11 
14 
20 
22 
27 
29 
30 #include "HLTDisplacedmumuFilter.h"
31 #include "TMath.h"
32 
33 //
34 // constructors and destructor
35 //
37  HLTFilter(iConfig),
38  fastAccept_ (iConfig.getParameter<bool>("FastAccept")),
39  minLxySignificance_ (iConfig.getParameter<double>("MinLxySignificance")),
40  maxLxySignificance_ (iConfig.getParameter<double>("MaxLxySignificance")),
41  maxNormalisedChi2_ (iConfig.getParameter<double>("MaxNormalisedChi2")),
42  minVtxProbability_ (iConfig.getParameter<double>("MinVtxProbability")),
43  minCosinePointingAngle_ (iConfig.getParameter<double>("MinCosinePointingAngle")),
44  DisplacedVertexTag_(iConfig.getParameter<edm::InputTag>("DisplacedVertexTag")),
45  DisplacedVertexToken_(consumes<reco::VertexCollection>(DisplacedVertexTag_)),
46  beamSpotTag_ (iConfig.getParameter<edm::InputTag> ("BeamSpotTag")),
47  beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_)),
48  MuonTag_ (iConfig.getParameter<edm::InputTag>("MuonTag")),
49  MuonToken_(consumes<reco::RecoChargedCandidateCollection>(MuonTag_))
50 {
51 }
52 
53 
55 
59  desc.add<bool>("FastAccept",false);
60  desc.add<double>("MinLxySignificance",0.0);
61  desc.add<double>("MaxLxySignificance",0.0);
62  desc.add<double>("MaxNormalisedChi2",10.0);
63  desc.add<double>("MinVtxProbability",0.0);
64  desc.add<double>("MinCosinePointingAngle",-2.0);
65  desc.add<edm::InputTag>("DisplacedVertexTag",edm::InputTag("hltDisplacedmumuVtx"));
66  desc.add<edm::InputTag>("BeamSpotTag",edm::InputTag("hltOnlineBeamSpot"));
67  desc.add<edm::InputTag>("MuonTag",edm::InputTag("hltL3MuonCandidates"));
68  descriptions.add("hltDisplacedmumuFilter", desc);
69  }
70 
71 // ------------ method called once each job just before starting event loop ------------
73 {
74 
75 }
76 
77 // ------------ method called once each job just after ending the event loop ------------
79 {
80 
81 }
82 
83 // ------------ method called on each new Event ------------
85 {
86 
87 
88  // get beam spot
89  reco::BeamSpot vertexBeamSpot;
90  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
91  iEvent.getByToken(beamSpotToken_,recoBeamSpotHandle);
92  vertexBeamSpot = *recoBeamSpotHandle;
93 
94 
95  // get displaced vertices
96  reco::VertexCollection displacedVertexColl;
97  edm::Handle<reco::VertexCollection> displacedVertexCollHandle;
98  bool foundVertexColl = iEvent.getByToken(DisplacedVertexToken_, displacedVertexCollHandle);
99  if(foundVertexColl) displacedVertexColl = *displacedVertexCollHandle;
100 
101 
102  // get muon collection
104  iEvent.getByToken(MuonToken_,mucands);
105 
106  // All HLT filters must create and fill an HLT filter object,
107  // recording any reconstructed physics objects satisfying (or not)
108  // this HLT filter, and place it in the Event.
109 
110 
111  // Ref to Candidate object to be recorded in filter object
114 
115  if (saveTags()) filterproduct.addCollectionTag(MuonTag_);
116 
117  bool triggered = false;
118 
119  // loop over vertex collection
120  for(auto displacedVertex : displacedVertexColl){
121  // check if the vertex actually consists of exactly two muon tracks, throw exception if not
122  if(displacedVertex.tracksSize() != 2) throw cms::Exception("BadLogic") << "HLTDisplacedmumuFilter: ERROR: the Jpsi vertex must have exactly two muons by definition. It now has n muons = "
123  << displacedVertex.tracksSize() << std::endl;
124 
125  float normChi2 = displacedVertex.normalizedChi2();
126  if (normChi2 > maxNormalisedChi2_) continue;
127 
128  double vtxProb = 0.0;
129  if( (displacedVertex.chi2()>=0.0) && (displacedVertex.ndof()>0) ) vtxProb = TMath::Prob(displacedVertex.chi2(), displacedVertex.ndof() );
130  if (vtxProb < minVtxProbability_) continue;
131 
132  // get the two muons from the vertex
133  auto trackIt = displacedVertex.tracks_begin();
134  reco::TrackRef vertextkRef1 = (*trackIt).castTo<reco::TrackRef>() ;
135  // the second one
136  trackIt++;
137  reco::TrackRef vertextkRef2 = (*trackIt).castTo<reco::TrackRef>();
138 
139  // first find these two tracks in the muon collection
140  reco::RecoChargedCandidateCollection::const_iterator cand1;
141  reco::RecoChargedCandidateCollection::const_iterator cand2;
142 
143  int iFoundRefs = 0;
144  for (auto cand=mucands->begin(); cand!=mucands->end(); cand++) {
145  reco::TrackRef tkRef = cand->get<reco::TrackRef>();
146  if(tkRef == vertextkRef1) {cand1 = cand; iFoundRefs++;}
147  if(tkRef == vertextkRef2) {cand2 = cand; iFoundRefs++;}
148  }
149  if(iFoundRefs < 2) throw cms::Exception("BadLogic") << "HLTDisplacedmumuFilter: ERROR: the muons matched with the Jpsi vertex tracks should be at least two by definition." << std::endl;
150 
151  // calculate two-track transverse momentum
152  math::XYZVector pperp(cand1->px() + cand2->px(),
153  cand1->py() + cand2->py(),
154  0.);
155 
156 
157  const reco::Vertex::Point& vpoint=displacedVertex.position();
158  //translate to global point, should be improved
159  GlobalPoint secondaryVertex (vpoint.x(), vpoint.y(), vpoint.z());
160 
161  reco::Vertex::Error verr = displacedVertex.error();
162  // translate to global error, should be improved
163  GlobalError err(verr.At(0,0), verr.At(1,0), verr.At(1,1), verr.At(2,0), verr.At(2,1), verr.At(2,2) );
164 
165  GlobalPoint displacementFromBeamspot( -1*((vertexBeamSpot.x0() - secondaryVertex.x()) + (secondaryVertex.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dxdz()),
166  -1*((vertexBeamSpot.y0() - secondaryVertex.y())+ (secondaryVertex.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dydz()), 0);
167 
168  float lxy = displacementFromBeamspot.perp();
169  float lxyerr = sqrt(err.rerr(displacementFromBeamspot));
170 
171 
172  //calculate the angle between the decay length and the mumu momentum
173  reco::Vertex::Point vperp(displacementFromBeamspot.x(),displacementFromBeamspot.y(),0.);
174 
175  float cosAlpha = vperp.Dot(pperp)/(vperp.R()*pperp.R());
176 
177  // check thresholds
178  if (cosAlpha < minCosinePointingAngle_) continue;
179  if (minLxySignificance_ > 0. && lxy/lxyerr < minLxySignificance_) continue;
180  if (maxLxySignificance_ > 0. && lxy/lxyerr > maxLxySignificance_) continue;
181  triggered = true;
182 
183  // now add the muons that passed to the filter object
184 
186  filterproduct.addObject(trigger::TriggerMuon,ref1);
187 
189  filterproduct.addObject(trigger::TriggerMuon,ref2);
190  }
191 
192  LogDebug("HLTDisplacedMumuFilter") << " >>>>> Result of HLTDisplacedMumuFilter is "<< triggered <<std::endl;
193 
194  return triggered;
195 }
196 
#define LogDebug(id)
HLTDisplacedmumuFilter(const edm::ParameterSet &)
double z0() const
z coordinate
Definition: BeamSpot.h:68
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
edm::Ref< RecoChargedCandidateCollection > RecoChargedCandidateRef
reference to an object in a collection of RecoChargedCandidate objects
math::Error< dimension >::type Error
covariance error matrix (3x3)
Definition: Vertex.h:43
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
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
~HLTDisplacedmumuFilter() override
int iEvent
Definition: GenABIO.cc:224
T sqrt(T t)
Definition: SSEVec.h:18
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
double dxdz() const
dxdz slope
Definition: BeamSpot.h:82
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
std::vector< RecoChargedCandidate > RecoChargedCandidateCollection
collectin of RecoChargedCandidate objects
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > MuonToken_
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
double y0() const
y coordinate
Definition: BeamSpot.h:66
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
edm::EDGetTokenT< reco::VertexCollection > DisplacedVertexToken_
double x0() const
x coordinate
Definition: BeamSpot.h:64