CMS 3D CMS Logo

HLTmumutkFilter.cc
Go to the documentation of this file.
1 #include <algorithm>
2 #include <cmath>
3 #include <vector>
4 
8 
15 
20 
21 #include "HLTmumutkFilter.h"
22 #include "TMath.h"
23 
24 using namespace edm;
25 using namespace reco;
26 using namespace std;
27 using namespace trigger;
28 
29 // ----------------------------------------------------------------------
31  muCandTag_ (iConfig.getParameter<edm::InputTag>("MuonTag")),
32  muCandToken_(consumes<reco::RecoChargedCandidateCollection>(muCandTag_)),
33  trkCandTag_ (iConfig.getParameter<edm::InputTag>("TrackTag")),
34  trkCandToken_(consumes<reco::RecoChargedCandidateCollection>(trkCandTag_)),
35  MuMuTkVertexTag_ (iConfig.getParameter<edm::InputTag>("MuMuTkVertexTag")),
36  MuMuTkVertexToken_(consumes<reco::VertexCollection>(MuMuTkVertexTag_)),
37  beamSpotTag_ (iConfig.getParameter<edm::InputTag> ("BeamSpotTag")),
38  beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_)),
39  maxEta_(iConfig.getParameter<double>("MaxEta")),
40  minPt_(iConfig.getParameter<double>("MinPt")),
41  maxNormalisedChi2_(iConfig.getParameter<double>("MaxNormalisedChi2")),
42  minVtxProbability_(iConfig.getParameter<double>("MinVtxProbability")),
43  minLxySignificance_(iConfig.getParameter<double>("MinLxySignificance")),
44  minCosinePointingAngle_(iConfig.getParameter<double>("MinCosinePointingAngle"))
45 {
46 }
47 
48 // ----------------------------------------------------------------------
50 
51 void
55  desc.add<double>("MaxEta",2.5);
56  desc.add<double>("MinPt" ,0.0);
57  desc.add<double>("MaxNormalisedChi2" ,10.0);
58  desc.add<double>("MinVtxProbability" , 0.0);
59  desc.add<double>("MinLxySignificance",3.0);
60  desc.add<double>("MinCosinePointingAngle",0.9);
61  desc.add<edm::InputTag>("MuonTag",edm::InputTag("hltL3MuonCandidates"));
62  desc.add<edm::InputTag>("TrackTag",edm::InputTag("hltMumukAllConeTracks"));
63  desc.add<edm::InputTag>("MuMuTkVertexTag",edm::InputTag("hltDisplacedmumuVtxProducerDoubleMu4Jpsi"));
64  desc.add<edm::InputTag>("BeamSpotTag",edm::InputTag("hltOnineBeamSpot"));
65  descriptions.add("HLTmumutkFilter",desc);
66 }
67 
68 // ----------------------------------------------------------------------
70 
71  //get the beamspot position
72  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
73  iEvent.getByToken(beamSpotToken_,recoBeamSpotHandle);
74  const reco::BeamSpot& vertexBeamSpot = *recoBeamSpotHandle;
75 
76  // get vertices
77  reco::VertexCollection displacedVertexColl;
78  edm::Handle<reco::VertexCollection> displacedVertexCollHandle;
79  bool foundVertexColl = iEvent.getByToken(MuMuTkVertexToken_, displacedVertexCollHandle);
80  if(foundVertexColl) displacedVertexColl = *displacedVertexCollHandle;
81 
82  // get muon collection
84  iEvent.getByToken(muCandToken_,mucands);
85 
86  // get track candidates around displaced muons
88  iEvent.getByToken(trkCandToken_,trkcands);
89 
90  // Ref to Candidate object to be recorded in filter object
94 
95  if (saveTags()) {
96  filterproduct.addCollectionTag(muCandTag_);
97  filterproduct.addCollectionTag(trkCandTag_);
98  }
99 
100  bool triggered = false;
101 
102  // loop over vertex collection
103  reco::VertexCollection::iterator it;
104  for(it = displacedVertexColl.begin(); it!= displacedVertexColl.end(); it++){
105  reco::Vertex displacedVertex = *it;
106 
107  // check if the vertex actually consists of exactly two muon + 1 track, throw exception if not
108  if(displacedVertex.tracksSize() != 3) throw cms::Exception("BadLogic") << "HLTmumutkFilter: ERROR: the Jpsi+trk vertex must have "
109  << "exactly two muons + 1 trk by definition. It now has n trakcs = "
110  << displacedVertex.tracksSize() << std::endl;
111 
112  float normChi2 = displacedVertex.normalizedChi2();
113  if (normChi2 > maxNormalisedChi2_) continue;
114 
115  double vtxProb = 0.0;
116  if ((displacedVertex.chi2()>=0.0) && (displacedVertex.ndof()>0) )
117  vtxProb = TMath::Prob(displacedVertex.chi2(), displacedVertex.ndof() );
118  if (vtxProb < minVtxProbability_) continue;
119 
120  // get the three tracks from the vertex
121  auto trackIt = displacedVertex.tracks_begin();
122  reco::TrackRef vertextkRef1 = (*trackIt).castTo<reco::TrackRef>() ;
123  trackIt++;
124  reco::TrackRef vertextkRef2 = (*trackIt).castTo<reco::TrackRef>();
125  trackIt++;
126  reco::TrackRef vertextkRef3 = (*trackIt).castTo<reco::TrackRef>();
127 
128  // first find the two muon tracks in the muon collection
129  reco::RecoChargedCandidateCollection::const_iterator mucand1;
130  reco::RecoChargedCandidateCollection::const_iterator mucand2;
131  reco::RecoChargedCandidateCollection::const_iterator tkcand ;
132 
133  int iFoundRefs = 0;
134  bool track1Matched = false;
135  bool track2Matched = false;
136  bool track3Matched = false;
137  for (auto cand=mucands->begin(); cand!=mucands->end(); cand++) {
138  reco::TrackRef tkRef = cand->get<reco::TrackRef>();
139  if (!track1Matched) {
140  if (tkRef == vertextkRef1 && iFoundRefs==0) {mucand1 = cand; iFoundRefs++; track1Matched = true;}
141  else if(tkRef == vertextkRef1 && iFoundRefs==1) {mucand2 = cand; iFoundRefs++; track1Matched = true;}
142  }
143  if (!track2Matched) {
144  if (tkRef == vertextkRef2 && iFoundRefs==0) {mucand1 = cand; iFoundRefs++; track2Matched = true;}
145  else if(tkRef == vertextkRef2 && iFoundRefs==1) {mucand2 = cand; iFoundRefs++; track2Matched = true;}
146  }
147  if (!track3Matched) {
148  if (tkRef == vertextkRef3 && iFoundRefs==0) {mucand1 = cand; iFoundRefs++; track3Matched = true;}
149  else if(tkRef == vertextkRef3 && iFoundRefs==1) {mucand2 = cand; iFoundRefs++; track3Matched = true;}
150  }
151  }
152  if(iFoundRefs < 2) throw cms::Exception("BadLogic") << "HLTmumutkFilterr: ERROR: the vertex must have "
153  << " at least two muons by definition." << std::endl;
154 
155  int iTrkFoundRefs = 0;
156  for (auto cand=trkcands->begin(); cand!=trkcands->end(); cand++) {
157  reco::TrackRef tkRef = cand->get<reco::TrackRef>();
158  if (tkRef == vertextkRef1 && iTrkFoundRefs==0 && !track1Matched) {tkcand = cand; iTrkFoundRefs++; track1Matched = true; break;}
159  if (tkRef == vertextkRef2 && iTrkFoundRefs==0 && !track2Matched) {tkcand = cand; iTrkFoundRefs++; track2Matched = true; break;}
160  if (tkRef == vertextkRef3 && iTrkFoundRefs==0 && !track3Matched) {tkcand = cand; iTrkFoundRefs++; track3Matched = true; break;}
161  }
162  if(iTrkFoundRefs == 0) throw cms::Exception("BadLogic") << "HLTmumutkFilterr: ERROR: the vertex must have "
163  << " at least one track by definition." << std::endl;
164 
165  // calculate three-track transverse momentum
166  math::XYZVector pperp(mucand1->px() + mucand2->px() + tkcand->px(),
167  mucand1->py() + mucand2->py() + tkcand->py(),
168  0.);
169 
170  // get vertex position and error to calculate the decay length significance
171  const reco::Vertex::Point& vpoint=displacedVertex.position();
172  reco::Vertex::Error verr = displacedVertex.error();
173  GlobalPoint secondaryVertex (vpoint.x(), vpoint.y(), vpoint.z());
174  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) );
175 
176  GlobalPoint displacementFromBeamspot( -1*((vertexBeamSpot.x0() - secondaryVertex.x()) + (secondaryVertex.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dxdz()),
177  -1*((vertexBeamSpot.y0() - secondaryVertex.y()) + (secondaryVertex.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dydz()),
178  0 );
179  float lxy = displacementFromBeamspot.perp();
180  float lxyerr = sqrt(err.rerr(displacementFromBeamspot));
181 
182  //calculate the angle between the decay length and the mumu momentum
183  Vertex::Point vperp(displacementFromBeamspot.x(),displacementFromBeamspot.y(),0.);
184  float cosAlpha = vperp.Dot(pperp)/(vperp.R()*pperp.R());
185 
186  if (pperp.R() < minPt_ ) continue;
187  if (lxy/lxyerr < minLxySignificance_ ) continue;
188  if (cosAlpha < minCosinePointingAngle_) continue;
189  triggered = true;
190 
191  refMu1=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (mucands,distance(mucands->begin(), mucand1)));
192  filterproduct.addObject(TriggerMuon,refMu1);
193  refMu2=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (mucands,distance(mucands->begin(), mucand2)));
194  filterproduct.addObject(TriggerMuon,refMu2);
195  refTrk=RecoChargedCandidateRef( Ref<RecoChargedCandidateCollection> (trkcands,distance(trkcands->begin(),tkcand)));
196  filterproduct.addObject(TriggerTrack,refTrk);
197 
198  }//end loop vertices
199 
200  LogDebug("HLTDisplacedMumuTrkFilter") << " >>>>> Result of HLTDisplacedMuMuTrkFilter is "<< triggered;
201  return triggered;
202 }
203 
204 
205 
206 bool HLTmumutkFilter::triggerdByPreviousLevel(const reco::RecoChargedCandidateRef & candref, const std::vector<reco::RecoChargedCandidateRef>& vcands){
207  unsigned int i=0;
208  unsigned int i_max=vcands.size();
209  for (;i!=i_max;++i){
210  if (candref == vcands[i]) return true;
211  }
212  return false;
213 }
#define LogDebug(id)
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
static bool triggerdByPreviousLevel(const reco::RecoChargedCandidateRef &, const std::vector< reco::RecoChargedCandidateRef > &)
double z0() const
z coordinate
Definition: BeamSpot.h:68
const double maxNormalisedChi2_
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
edm::InputTag muCandTag_
HLTmumutkFilter(const edm::ParameterSet &)
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
const double minLxySignificance_
const Point & position() const
position
Definition: Vertex.h:109
edm::InputTag trkCandTag_
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
~HLTmumutkFilter() override
double dydz() const
dydz slope
Definition: BeamSpot.h:84
int iEvent
Definition: GenABIO.cc:224
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > trkCandToken_
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > muCandToken_
T sqrt(T t)
Definition: SSEVec.h:18
const double minCosinePointingAngle_
double chi2() const
chi-squares
Definition: Vertex.h:98
const double minVtxProbability_
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
ParameterDescriptionBase * add(U const &iLabel, T const &value)
const double minPt_
double dxdz() const
dxdz slope
Definition: BeamSpot.h:82
double ndof() const
Definition: Vertex.h:105
edm::EDGetTokenT< reco::VertexCollection > MuMuTkVertexToken_
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
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
Error error() const
return SMatrix
Definition: Vertex.h:139
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
double y0() const
y coordinate
Definition: BeamSpot.h:66
double normalizedChi2() const
chi-squared divided by n.d.o.f.
Definition: Vertex.h:107
trackRef_iterator tracks_begin() const
first iterator over tracks
Definition: Vertex.cc:76
size_t tracksSize() const
number of tracks
Definition: Vertex.cc:71
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
double x0() const
x coordinate
Definition: BeamSpot.h:64