btau/plugins/HLTDisplacedtktkVtxProducer.cc

0001 #include <iostream>
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0005 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0006 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0007 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0008 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0009 #include "DataFormats/Candidate/interface/Candidate.h"
0010 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0011 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0012 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0013 #include "DataFormats/TrackReco/interface/Track.h"
0014 #include "DataFormats/VertexReco/interface/Vertex.h"
0015 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0016 #include "DataFormats/Common/interface/RefToBase.h"
0017 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
0018 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0019 #include "HLTDisplacedtktkVtxProducer.h"
0020
0021 using namespace edm;
0022 using namespace reco;
0023 using namespace std;
0024 using namespace trigger;
0025 //
0026 // constructors and destructor
0027 //
0028 HLTDisplacedtktkVtxProducer::HLTDisplacedtktkVtxProducer(const edm::ParameterSet& iConfig)
0029     : transientTrackRecordToken_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))),
0030       srcTag_(iConfig.getParameter<edm::InputTag>("Src")),
0031       srcToken_(consumes<reco::RecoChargedCandidateCollection>(srcTag_)),
0032       previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
0033       previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0034       maxEta_(iConfig.getParameter<double>("MaxEta")),
0035       minPt_(iConfig.getParameter<double>("MinPt")),
0036       minPtPair_(iConfig.getParameter<double>("MinPtPair")),
0037       minInvMass_(iConfig.getParameter<double>("MinInvMass")),
0038       maxInvMass_(iConfig.getParameter<double>("MaxInvMass")),
0039       massParticle1_(iConfig.getParameter<double>("massParticle1")),
0040       massParticle2_(iConfig.getParameter<double>("massParticle2")),
0041       chargeOpt_(iConfig.getParameter<int>("ChargeOpt")),
0042       triggerTypeDaughters_(iConfig.getParameter<int>("triggerTypeDaughters")) {
0043   produces<VertexCollection>();
0044 }
0045
0046 HLTDisplacedtktkVtxProducer::~HLTDisplacedtktkVtxProducer() = default;
0047
0048 void HLTDisplacedtktkVtxProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0049   edm::ParameterSetDescription desc;
0050   desc.add<edm::InputTag>("Src", edm::InputTag("hltL3MuonCandidates"));
0051   desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
0052   desc.add<double>("MaxEta", 2.5);
0053   desc.add<double>("MinPt", 0.0);
0054   desc.add<double>("MinPtPair", 0.0);
0055   desc.add<double>("MinInvMass", 1.0);
0056   desc.add<double>("MaxInvMass", 20.0);
0057   desc.add<double>("massParticle1", 0.1396);
0058   desc.add<double>("massParticle2", 0.4937);
0059   desc.add<int>("ChargeOpt", -1);
0060   desc.add<int>("triggerTypeDaughters", 0);
0061
0062   descriptions.add("hltDisplacedtktkVtxProducer", desc);
0063 }
0064
0065 // ------------ method called on each new Event  ------------
0066 void HLTDisplacedtktkVtxProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0067   double const firstTrackMass = massParticle1_;
0068   double const firstTrackMass2 = firstTrackMass * firstTrackMass;
0069   double const secondTrackMass = massParticle2_;
0070   double const secondTrackMass2 = secondTrackMass * secondTrackMass;
0071
0072   // get hold of track trks
0073   Handle<RecoChargedCandidateCollection> trackcands;
0074   iEvent.getByToken(srcToken_, trackcands);
0075
0076   // get the transient track builder
0077   auto const& theB = iSetup.getHandle(transientTrackRecordToken_);
0078
0079   std::unique_ptr<VertexCollection> vertexCollection(new VertexCollection());
0080
0081   // look at all trackcands,  check cuts and make vertices
0082   double e1, e2;
0083   Particle::LorentzVector p, p1, p2;
0084
0085   RecoChargedCandidateCollection::const_iterator cand1;
0086   RecoChargedCandidateCollection::const_iterator cand2;
0087
0088   // get the objects passing the previous filter
0089   Handle<TriggerFilterObjectWithRefs> previousCands;
0090   iEvent.getByToken(previousCandToken_, previousCands);
0091
0092   vector<RecoChargedCandidateRef> vPrevCands;
0093   previousCands->getObjects(triggerTypeDaughters_, vPrevCands);
0094
0095   std::vector<bool> candComp;
0096   for (cand1 = trackcands->begin(); cand1 != trackcands->end(); cand1++)
0097     candComp.push_back(checkPreviousCand(cand1->get<TrackRef>(), vPrevCands));
0098
0099   for (cand1 = trackcands->begin(); cand1 != trackcands->end(); cand1++) {
0100     TrackRef tk1 = cand1->get<TrackRef>();
0101     LogDebug("HLTDisplacedtktkVtxProducer") << " 1st track in loop: q*pt= " << cand1->charge() * cand1->pt()
0102                                             << ", eta= " << cand1->eta() << ", hits= " << tk1->numberOfValidHits();
0103
0104     //first check if this track passed the previous filter
0105     if (!candComp[cand1 - trackcands->begin()])
0106       continue;
0107     // if( ! checkPreviousCand( tk1, vPrevCands) ) continue;
0108
0109     // cuts
0110     if (abs(cand1->eta()) > maxEta_)
0111       continue;
0112     if (cand1->pt() < minPt_)
0113       continue;
0114
0115     cand2 = trackcands->begin();
0116     if (massParticle1_ == massParticle2_) {
0117       cand2 = cand1 + 1;
0118     }
0119
0120     for (; cand2 != trackcands->end(); cand2++) {
0121       TrackRef tk2 = cand2->get<TrackRef>();
0122       if (tk1 == tk2)
0123         continue;
0124
0125       // eta cut
0126       LogDebug("HLTDisplacedtktkVtxProducer")
0127           << " 2nd track in loop: q*pt= " << cand2->charge() * cand2->pt() << ", eta= " << cand2->eta()
0128           << ", hits= " << tk2->numberOfValidHits() << ", d0= " << tk2->d0();
0129       //first check if this track passed the previous filter
0130       if (!candComp[cand2 - trackcands->begin()])
0131         continue;
0132       // if( ! checkPreviousCand( tk2, vPrevCands) ) continue;
0133
0134       // cuts
0135       if (abs(cand2->eta()) > maxEta_)
0136         continue;
0137       if (cand2->pt() < minPt_)
0138         continue;
0139
0140       // opposite sign or same sign
0141       if (chargeOpt_ < 0) {
0142         if (cand1->charge() * cand2->charge() > 0)
0143           continue;
0144       } else if (chargeOpt_ > 0) {
0145         if (cand1->charge() * cand2->charge() < 0)
0146           continue;
0147       }
0148
0149       // Combined ditrack system
0150       e1 = sqrt(cand1->momentum().Mag2() + firstTrackMass2);
0151       e2 = sqrt(cand2->momentum().Mag2() + secondTrackMass2);
0152       p1 = Particle::LorentzVector(cand1->px(), cand1->py(), cand1->pz(), e1);
0153       p2 = Particle::LorentzVector(cand2->px(), cand2->py(), cand2->pz(), e2);
0154       p = p1 + p2;
0155
0156       if (p.pt() < minPtPair_)
0157         continue;
0158
0159       double invmass = abs(p.mass());
0160       LogDebug("HLTDisplacedtktkVtxProducer") << " ... 1-2 invmass= " << invmass;
0161
0162       if (invmass < minInvMass_)
0163         continue;
0164       if (invmass > maxInvMass_)
0165         continue;
0166
0167       // do the vertex fit
0168       vector<TransientTrack> t_tks;
0169       TransientTrack ttkp1 = (*theB).build(&tk1);
0170       TransientTrack ttkp2 = (*theB).build(&tk2);
0171       t_tks.push_back(ttkp1);
0172       t_tks.push_back(ttkp2);
0173
0174       if (t_tks.size() != 2)
0175         continue;
0176
0177       KalmanVertexFitter kvf;
0178       TransientVertex tv = kvf.vertex(t_tks);
0179
0180       if (!tv.isValid())
0181         continue;
0182
0183       Vertex vertex = tv;
0184
0185       // put vertex in the event
0186       vertexCollection->push_back(vertex);
0187     }
0188   }
0189   iEvent.put(std::move(vertexCollection));
0190 }
0191
0192 bool HLTDisplacedtktkVtxProducer::checkPreviousCand(const TrackRef& trackref,
0193                                                     const vector<RecoChargedCandidateRef>& refVect) const {
0194   bool ok = false;
0195   for (auto& i : refVect) {
0196     if (i->get<TrackRef>() == trackref) {
0197       ok = true;
0198       break;
0199     }
0200   }
0201   return ok;
0202 }