CMS 3D CMS Logo

HLTDisplacedtktktkFilter.cc
Go to the documentation of this file.
1 #include <iostream>
2 
5 
11 
14 
20 
22 
27 
29 
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  triggerTypeDaughters_(iConfig.getParameter<int>("triggerTypeDaughters")),
45  DisplacedVertexTag_(iConfig.getParameter<edm::InputTag>("DisplacedVertexTag")),
46  DisplacedVertexToken_(consumes<reco::VertexCollection>(DisplacedVertexTag_)),
47  beamSpotTag_ (iConfig.getParameter<edm::InputTag> ("BeamSpotTag")),
48  beamSpotToken_(consumes<reco::BeamSpot>(beamSpotTag_)),
49  TrackTag_ (iConfig.getParameter<edm::InputTag>("TrackTag")),
50  TrackToken_(consumes<reco::RecoChargedCandidateCollection>(TrackTag_))
51 {
52 }
53 
54 
56 
60  desc.add<bool>("FastAccept",false);
61  desc.add<double>("MinLxySignificance",0.0);
62  desc.add<double>("MaxLxySignificance",0.0);
63  desc.add<double>("MaxNormalisedChi2",10.0);
64  desc.add<double>("MinVtxProbability",0.0);
65  desc.add<double>("MinCosinePointingAngle",-2.0);
66  desc.add<int>("triggerTypeDaughters",0);
67  desc.add<edm::InputTag>("DisplacedVertexTag",edm::InputTag("hltDisplacedtktktkVtx"));
68  desc.add<edm::InputTag>("BeamSpotTag",edm::InputTag("hltOnlineBeamSpot"));
69  desc.add<edm::InputTag>("TrackTag",edm::InputTag("hltL3MuonCandidates"));
70  descriptions.add("hltDisplacedtktktkFilter", desc);
71 }
72 
73 // ------------ method called once each job just before starting event loop ------------
75 {
76 
77 }
78 
79 // ------------ method called once each job just after ending the event loop ------------
81 {
82 
83 }
84 
85 // ------------ method called on each new Event ------------
87 {
88 
89 
90  // get beam spot
91  reco::BeamSpot vertexBeamSpot;
92  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
93  iEvent.getByToken(beamSpotToken_,recoBeamSpotHandle);
94  vertexBeamSpot = *recoBeamSpotHandle;
95 
96 
97  // get displaced vertices
98  reco::VertexCollection displacedVertexColl;
99  edm::Handle<reco::VertexCollection> displacedVertexCollHandle;
100  bool foundVertexColl = iEvent.getByToken(DisplacedVertexToken_, displacedVertexCollHandle);
101  if(foundVertexColl) displacedVertexColl = *displacedVertexCollHandle;
102 
103 
104  // get track collection
106  iEvent.getByToken(TrackToken_,trackcands);
107 
108  // All HLT filters must create and fill an HLT filter object,
109  // recording any reconstructed physics objects satisfying (or not)
110  // this HLT filter, and place it in the Event.
111 
112 
113  if (saveTags()) filterproduct.addCollectionTag(TrackTag_);
114 
115  bool triggered = false;
116 
117  // loop over vertex collection
118  for(auto const &displacedVertex : displacedVertexColl){
119  // check if the vertex actually consists of exactly three track tracks, reject the event if not
120  if(displacedVertex.tracksSize() != 3){
121  edm::LogError("HLTDisplacedtktktkFilter") << "HLTDisplacedtktktkFilter: ERROR: the tktktk vertex must have exactly three tracks by definition. It now has n tracks = "<< displacedVertex.tracksSize() << std::endl;
122  return false;
123  }
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 three tracks 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  // the third one
139  trackIt++;
140  reco::TrackRef vertextkRef3 = (*trackIt).castTo<reco::TrackRef>();
141 
142  // first find these three tracks in the track collection
143  reco::RecoChargedCandidateCollection::const_iterator cand1;
144  reco::RecoChargedCandidateCollection::const_iterator cand2;
145  reco::RecoChargedCandidateCollection::const_iterator cand3;
146 
147  int iFoundRefs = 0;
148  for (auto cand=trackcands->begin(); cand!=trackcands->end(); cand++) {
149  reco::TrackRef tkRef = cand->get<reco::TrackRef>();
150  if(tkRef == vertextkRef1) {cand1 = cand; iFoundRefs++;}
151  if(tkRef == vertextkRef2) {cand2 = cand; iFoundRefs++;}
152  if(tkRef == vertextkRef3) {cand3 = cand; iFoundRefs++;}
153  }
154 
155  if(iFoundRefs < 3){
156  edm::LogError("HLTDisplacedtktktkFilter") << "HLTDisplacedtktktkFilter: ERROR: the tracks matched with the tktktk vertex tracks should be at least three by definition." << std::endl;
157  return false;
158  }
159  // calculate three-track transverse momentum
160  math::XYZVector pperp(cand1->px() + cand2->px() + cand3->px(),
161  cand1->py() + cand2->py() + cand3->py(),
162  0.);
163 
164  const reco::Vertex::Point& vpoint=displacedVertex.position();
165  reco::Vertex::Error verr = displacedVertex.error();
166  // translate to global error, should be improved
167  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) );
168 
169  GlobalPoint displacementFromBeamspot( -1*((vertexBeamSpot.x0() - vpoint.x()) + (vpoint.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dxdz()),
170  -1*((vertexBeamSpot.y0() - vpoint.y())+ (vpoint.z() - vertexBeamSpot.z0()) * vertexBeamSpot.dydz()), 0);
171 
172  float lxy = displacementFromBeamspot.perp();
173  float lxyerr = sqrt(err.rerr(displacementFromBeamspot));
174 
175 
176  //calculate the angle between the decay length and the tktk momentum
177  reco::Vertex::Point vperp(displacementFromBeamspot.x(),displacementFromBeamspot.y(),0.);
178 
179  float cosAlpha = vperp.Dot(pperp)/(vperp.R()*pperp.R());
180 
181  // check thresholds
182  if (cosAlpha < minCosinePointingAngle_) continue;
183  if (minLxySignificance_ > 0. && lxy/lxyerr < minLxySignificance_) continue;
184  if (maxLxySignificance_ > 0. && lxy/lxyerr > maxLxySignificance_) continue;
185  triggered = true;
186 
187  // now add the tracks that passed to the filter object
188 
189  // Ref to Candidate object to be recorded in filter object
193 
194  ref1=reco::RecoChargedCandidateRef( edm::Ref<reco::RecoChargedCandidateCollection> (trackcands,distance(trackcands->begin(), cand1)));
195  filterproduct.addObject(triggerTypeDaughters_,ref1);
196 
197  ref2=reco::RecoChargedCandidateRef( edm::Ref<reco::RecoChargedCandidateCollection> (trackcands,distance(trackcands->begin(),cand2 )));
198  filterproduct.addObject(triggerTypeDaughters_,ref2);
199 
200  ref3=reco::RecoChargedCandidateRef( edm::Ref<reco::RecoChargedCandidateCollection> (trackcands,distance(trackcands->begin(),cand3 )));
201  filterproduct.addObject(triggerTypeDaughters_,ref3);
202  }
203 
204  LogDebug("HLTDisplacedtktktkFilter") << " >>>>> Result of HLTDisplacedtktktkFilter is "<< triggered <<std::endl;
205 
206  return triggered;
207 }
#define LogDebug(id)
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
~HLTDisplacedtktktkFilter() override
edm::EDGetTokenT< reco::BeamSpot > beamSpotToken_
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
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
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
HLTDisplacedtktktkFilter(const edm::ParameterSet &)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
double dxdz() const
dxdz slope
Definition: BeamSpot.h:82
std::vector< RecoChargedCandidate > RecoChargedCandidateCollection
collectin of RecoChargedCandidate objects
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > TrackToken_
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)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
double y0() const
y coordinate
Definition: BeamSpot.h:66
edm::EDGetTokenT< reco::VertexCollection > DisplacedVertexToken_
double x0() const
x coordinate
Definition: BeamSpot.h:64
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)