File indexing completed on 2024-04-06 12:15:49
0001 #include <iostream>
0002
0003 #include "FWCore/Framework/interface/ESHandle.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0006 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0007 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0008 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0009 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0010 #include "DataFormats/Candidate/interface/Candidate.h"
0011 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0012 #include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"
0013 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0014 #include "DataFormats/TrackReco/interface/Track.h"
0015 #include "DataFormats/VertexReco/interface/Vertex.h"
0016 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0017 #include "DataFormats/Common/interface/RefToBase.h"
0018 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
0019 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0020 #include "HLTDisplacedmumuVtxProducer.h"
0021
0022 using namespace edm;
0023 using namespace reco;
0024 using namespace std;
0025 using namespace trigger;
0026
0027
0028
0029 HLTDisplacedmumuVtxProducer::HLTDisplacedmumuVtxProducer(const edm::ParameterSet& iConfig)
0030 : transientTrackRecordToken_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))),
0031 srcTag_(iConfig.getParameter<edm::InputTag>("Src")),
0032 srcToken_(consumes<reco::RecoChargedCandidateCollection>(srcTag_)),
0033 previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
0034 previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
0035 matchToPrevious_(iConfig.getParameter<bool>("matchToPrevious")),
0036 maxEta_(iConfig.getParameter<double>("MaxEta")),
0037 minPt_(iConfig.getParameter<double>("MinPt")),
0038 minPtPair_(iConfig.getParameter<double>("MinPtPair")),
0039 minInvMass_(iConfig.getParameter<double>("MinInvMass")),
0040 maxInvMass_(iConfig.getParameter<double>("MaxInvMass")),
0041 chargeOpt_(iConfig.getParameter<int>("ChargeOpt")) {
0042 produces<VertexCollection>();
0043 }
0044
0045 HLTDisplacedmumuVtxProducer::~HLTDisplacedmumuVtxProducer() = default;
0046
0047 void HLTDisplacedmumuVtxProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0048 edm::ParameterSetDescription desc;
0049 desc.add<edm::InputTag>("Src", edm::InputTag("hltL3MuonCandidates"));
0050 desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
0051 desc.add<bool>("matchToPrevious", true);
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<int>("ChargeOpt", -1);
0058 descriptions.add("hltDisplacedmumuVtxProducer", desc);
0059 }
0060
0061
0062 void HLTDisplacedmumuVtxProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0063 double const MuMass = 0.106;
0064 double const MuMass2 = MuMass * MuMass;
0065
0066
0067 Handle<RecoChargedCandidateCollection> mucands;
0068 iEvent.getByToken(srcToken_, mucands);
0069
0070
0071 auto const& theB = iSetup.getHandle(transientTrackRecordToken_);
0072
0073 std::unique_ptr<VertexCollection> vertexCollection(new VertexCollection());
0074
0075
0076 double e1, e2;
0077 Particle::LorentzVector p, p1, p2;
0078
0079 RecoChargedCandidateCollection::const_iterator cand1;
0080 RecoChargedCandidateCollection::const_iterator cand2;
0081
0082
0083 Handle<TriggerFilterObjectWithRefs> previousCands;
0084 if (matchToPrevious_)
0085 iEvent.getByToken(previousCandToken_, previousCands);
0086
0087 vector<RecoChargedCandidateRef> vPrevCands;
0088 if (matchToPrevious_)
0089 previousCands->getObjects(TriggerMuon, vPrevCands);
0090
0091 for (cand1 = mucands->begin(); cand1 != mucands->end(); cand1++) {
0092 TrackRef tk1 = cand1->get<TrackRef>();
0093 LogDebug("HLTDisplacedMumuFilter") << " 1st muon in loop: q*pt= " << cand1->charge() * cand1->pt()
0094 << ", eta= " << cand1->eta() << ", hits= " << tk1->numberOfValidHits();
0095
0096
0097 if (matchToPrevious_ && !checkPreviousCand(tk1, vPrevCands))
0098 continue;
0099
0100
0101 if (fabs(cand1->eta()) > maxEta_)
0102 continue;
0103 if (cand1->pt() < minPt_)
0104 continue;
0105
0106 cand2 = cand1;
0107 cand2++;
0108 for (; cand2 != mucands->end(); cand2++) {
0109 TrackRef tk2 = cand2->get<TrackRef>();
0110
0111
0112 LogDebug("HLTDisplacedmumuVtxProducer")
0113 << " 2nd muon in loop: q*pt= " << cand2->charge() * cand2->pt() << ", eta= " << cand2->eta()
0114 << ", hits= " << tk2->numberOfValidHits() << ", d0= " << tk2->d0();
0115
0116 if (matchToPrevious_ && !checkPreviousCand(tk2, vPrevCands))
0117 continue;
0118
0119
0120 if (fabs(cand2->eta()) > maxEta_)
0121 continue;
0122 if (cand2->pt() < minPt_)
0123 continue;
0124
0125
0126 if (chargeOpt_ < 0) {
0127 if (cand1->charge() * cand2->charge() > 0)
0128 continue;
0129 } else if (chargeOpt_ > 0) {
0130 if (cand1->charge() * cand2->charge() < 0)
0131 continue;
0132 }
0133
0134
0135 e1 = sqrt(cand1->momentum().Mag2() + MuMass2);
0136 e2 = sqrt(cand2->momentum().Mag2() + MuMass2);
0137 p1 = Particle::LorentzVector(cand1->px(), cand1->py(), cand1->pz(), e1);
0138 p2 = Particle::LorentzVector(cand2->px(), cand2->py(), cand2->pz(), e2);
0139 p = p1 + p2;
0140
0141 if (p.pt() < minPtPair_)
0142 continue;
0143
0144 double invmass = abs(p.mass());
0145 LogDebug("HLTDisplacedMumuFilter") << " ... 1-2 invmass= " << invmass;
0146
0147 if (invmass < minInvMass_)
0148 continue;
0149 if (invmass > maxInvMass_)
0150 continue;
0151
0152
0153 vector<TransientTrack> t_tks;
0154 TransientTrack ttkp1 = (*theB).build(&tk1);
0155 TransientTrack ttkp2 = (*theB).build(&tk2);
0156 t_tks.push_back(ttkp1);
0157 t_tks.push_back(ttkp2);
0158
0159 if (t_tks.size() != 2)
0160 continue;
0161
0162 KalmanVertexFitter kvf;
0163 TransientVertex tv = kvf.vertex(t_tks);
0164
0165 if (!tv.isValid())
0166 continue;
0167
0168 Vertex vertex = tv;
0169
0170
0171 vertexCollection->push_back(vertex);
0172 }
0173 }
0174 iEvent.put(std::move(vertexCollection));
0175 }
0176
0177 bool HLTDisplacedmumuVtxProducer::checkPreviousCand(const TrackRef& trackref,
0178 const vector<RecoChargedCandidateRef>& refVect) const {
0179 bool ok = false;
0180 for (auto& i : refVect) {
0181 if (i->get<TrackRef>() == trackref) {
0182 ok = true;
0183 break;
0184 }
0185 }
0186 return ok;
0187 }