Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:49:44

0001 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0002 #include "HLTDoubletDZ.h"
0003 
0004 #include "DataFormats/Common/interface/Handle.h"
0005 
0006 #include "DataFormats/Common/interface/Ref.h"
0007 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0008 
0009 //#include "DataFormats/Candidate/interface/Particle.h"
0010 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0011 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateFwd.h"
0012 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0013 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
0014 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0015 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0016 #include "DataFormats/L1TCorrelator/interface/TkMuon.h"
0017 #include "DataFormats/L1TCorrelator/interface/TkMuonFwd.h"
0018 #include "DataFormats/L1TParticleFlow/interface/PFTau.h"
0019 #include "DataFormats/L1TParticleFlow/interface/HPSPFTau.h"
0020 #include "DataFormats/L1TParticleFlow/interface/HPSPFTauFwd.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 
0023 #include "FWCore/Framework/interface/MakerMacros.h"
0024 
0025 #include "DataFormats/Math/interface/deltaR.h"
0026 
0027 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0028 #include <cmath>
0029 
0030 //
0031 // constructors and destructor
0032 //
0033 template <typename T1, typename T2>
0034 HLTDoubletDZ<T1, T2>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0035     : HLTFilter(iConfig),
0036       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0037       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0038       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0039       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0040       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0041       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0042       //electronToken_ (consumes<reco::ElectronCollection>(iConfig.template getParameter<edm::InputTag>("electronTag"))),
0043       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0044       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0045       minDR_(iConfig.template getParameter<double>("MinDR")),
0046       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0047       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0048       min_N_(iConfig.template getParameter<int>("MinN")),
0049       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0050       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0051 {}
0052 
0053 template <>
0054 HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoEcalCandidate>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0055     : HLTFilter(iConfig),
0056       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0057       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0058       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0059       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0060       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0061       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0062       electronToken_(consumes<reco::ElectronCollection>(iConfig.template getParameter<edm::InputTag>("electronTag"))),
0063       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0064       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0065       minDR_(iConfig.template getParameter<double>("MinDR")),
0066       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0067       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0068       min_N_(iConfig.template getParameter<int>("MinN")),
0069       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0070       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0071 {}
0072 
0073 template <>
0074 HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoEcalCandidate>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0075     : HLTFilter(iConfig),
0076       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0077       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0078       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0079       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0080       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0081       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0082       electronToken_(consumes<reco::ElectronCollection>(iConfig.template getParameter<edm::InputTag>("electronTag"))),
0083       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0084       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0085       minDR_(iConfig.template getParameter<double>("MinDR")),
0086       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0087       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0088       min_N_(iConfig.template getParameter<int>("MinN")),
0089       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0090       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0091 {}
0092 
0093 template <>
0094 HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoChargedCandidate>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0095     : HLTFilter(iConfig),
0096       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0097       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0098       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0099       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0100       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0101       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0102       electronToken_(consumes<reco::ElectronCollection>(iConfig.template getParameter<edm::InputTag>("electronTag"))),
0103       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0104       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0105       minDR_(iConfig.template getParameter<double>("MinDR")),
0106       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0107       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0108       min_N_(iConfig.template getParameter<int>("MinN")),
0109       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0110       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0111 {}
0112 
0113 template <>
0114 HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoChargedCandidate>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0115     : HLTFilter(iConfig),
0116       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0117       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0118       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0119       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0120       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0121       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0122       //electronToken_ (consumes<reco::ElectronCollection>(iConfig.template getParameter<edm::InputTag>("electronTag"))),
0123       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0124       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0125       minDR_(iConfig.template getParameter<double>("MinDR")),
0126       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0127       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0128       min_N_(iConfig.template getParameter<int>("MinN")),
0129       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0130       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0131 {}
0132 
0133 template <>
0134 HLTDoubletDZ<l1t::TkMuon, l1t::TkMuon>::HLTDoubletDZ(const edm::ParameterSet& iConfig)
0135     : HLTFilter(iConfig),
0136       originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0137       originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0138       inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0139       inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0140       inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0141       inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0142       triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0143       triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0144       minDR_(iConfig.template getParameter<double>("MinDR")),
0145       maxDZ_(iConfig.template getParameter<double>("MaxDZ")),
0146       minPixHitsForDZ_(iConfig.template getParameter<int>("MinPixHitsForDZ")),
0147       min_N_(iConfig.template getParameter<int>("MinN")),
0148       checkSC_(iConfig.template getParameter<bool>("checkSC")),
0149       same_(inputTag1_.encode() == inputTag2_.encode())  // same collections to be compared?
0150 {}
0151 
0152 template <typename T1, typename T2>
0153 HLTDoubletDZ<T1, T2>::~HLTDoubletDZ() {}
0154 
0155 template <typename T1, typename T2>
0156 void HLTDoubletDZ<T1, T2>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0157   edm::ParameterSetDescription desc;
0158   makeHLTFilterDescription(desc);
0159   std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0160   std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0161   desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0162   desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0163   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0164   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0165   desc.add<int>("triggerType1", 0);
0166   desc.add<int>("triggerType2", 0);
0167   desc.add<double>("MinDR", -1.0);
0168   desc.add<double>("MaxDZ", 0.2);
0169   desc.add<int>("MinPixHitsForDZ", 0);
0170   desc.add<bool>("checkSC", false);
0171   desc.add<int>("MinN", 1);
0172   descriptions.add(defaultModuleLabel<HLTDoubletDZ<T1, T2>>(), desc);
0173 }
0174 
0175 template <>
0176 void HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoEcalCandidate>::fillDescriptions(
0177     edm::ConfigurationDescriptions& descriptions) {
0178   edm::ParameterSetDescription desc;
0179   makeHLTFilterDescription(desc);
0180   std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0181   std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0182   desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0183   desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0184   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0185   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0186   desc.add<edm::InputTag>("electronTag", edm::InputTag("electronTag"));
0187   desc.add<int>("triggerType1", 0);
0188   desc.add<int>("triggerType2", 0);
0189   desc.add<double>("MinDR", -1.0);
0190   desc.add<double>("MaxDZ", 0.2);
0191   desc.add<int>("MinPixHitsForDZ", 0);
0192   desc.add<bool>("checkSC", false);
0193   desc.add<int>("MinN", 1);
0194   descriptions.add(defaultModuleLabel<HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoEcalCandidate>>(), desc);
0195 }
0196 
0197 template <>
0198 void HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoEcalCandidate>::fillDescriptions(
0199     edm::ConfigurationDescriptions& descriptions) {
0200   edm::ParameterSetDescription desc;
0201   makeHLTFilterDescription(desc);
0202   std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0203   std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0204   desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0205   desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0206   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0207   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0208   desc.add<edm::InputTag>("electronTag", edm::InputTag("electronTag"));
0209   desc.add<int>("triggerType1", 0);
0210   desc.add<int>("triggerType2", 0);
0211   desc.add<double>("MinDR", -1.0);
0212   desc.add<double>("MaxDZ", 0.2);
0213   desc.add<int>("MinPixHitsForDZ", 0);
0214   desc.add<bool>("checkSC", false);
0215   desc.add<int>("MinN", 1);
0216   descriptions.add(defaultModuleLabel<HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoEcalCandidate>>(), desc);
0217 }
0218 
0219 template <>
0220 void HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoChargedCandidate>::fillDescriptions(
0221     edm::ConfigurationDescriptions& descriptions) {
0222   edm::ParameterSetDescription desc;
0223   makeHLTFilterDescription(desc);
0224   std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0225   std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0226   desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0227   desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0228   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0229   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0230   desc.add<edm::InputTag>("electronTag", edm::InputTag("electronTag"));
0231   desc.add<int>("triggerType1", 0);
0232   desc.add<int>("triggerType2", 0);
0233   desc.add<double>("MinDR", -1.0);
0234   desc.add<double>("MaxDZ", 0.2);
0235   desc.add<int>("MinPixHitsForDZ", 0);
0236   desc.add<bool>("checkSC", false);
0237   desc.add<int>("MinN", 1);
0238   descriptions.add(defaultModuleLabel<HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoChargedCandidate>>(), desc);
0239 }
0240 
0241 template <>
0242 void HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoChargedCandidate>::fillDescriptions(
0243     edm::ConfigurationDescriptions& descriptions) {
0244   edm::ParameterSetDescription desc;
0245   makeHLTFilterDescription(desc);
0246   std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0247   std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0248   desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0249   desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0250   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0251   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0252   desc.add<int>("triggerType1", 0);
0253   desc.add<int>("triggerType2", 0);
0254   desc.add<double>("MinDR", -1.0);
0255   desc.add<double>("MaxDZ", 0.2);
0256   desc.add<int>("MinPixHitsForDZ", 0);
0257   desc.add<bool>("checkSC", false);
0258   desc.add<int>("MinN", 1);
0259   descriptions.add(defaultModuleLabel<HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoChargedCandidate>>(), desc);
0260 }
0261 
0262 template <typename T1, typename T2>
0263 bool
0264 //HLTDoubletDZ<T1, T2>::getCollections(edm::Event& iEvent, std::vector<T1Ref>& coll1, std::vector<T2Ref>& coll2) const {
0265 HLTDoubletDZ<T1, T2>::getCollections(edm::Event& iEvent,
0266                                      std::vector<T1Ref>& coll1,
0267                                      std::vector<T2Ref>& coll2,
0268                                      trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0269   edm::Handle<trigger::TriggerFilterObjectWithRefs> handle1, handle2;
0270   if (iEvent.getByToken(inputToken1_, handle1) and iEvent.getByToken(inputToken2_, handle2)) {
0271     // get hold of pre-filtered object collections
0272     handle1->getObjects(triggerType1_, coll1);
0273     handle2->getObjects(triggerType2_, coll2);
0274     const trigger::size_type n1(coll1.size());
0275     const trigger::size_type n2(coll2.size());
0276 
0277     if (saveTags()) {
0278       edm::InputTag tagOld;
0279       for (unsigned int i = 0; i < originTag1_.size(); ++i) {
0280         filterproduct.addCollectionTag(originTag1_[i]);
0281       }
0282       tagOld = edm::InputTag();
0283       for (trigger::size_type i1 = 0; i1 != n1; ++i1) {
0284         const edm::ProductID pid(coll1[i1].id());
0285         const auto& prov = iEvent.getStableProvenance(pid);
0286         const std::string& label(prov.moduleLabel());
0287         const std::string& instance(prov.productInstanceName());
0288         const std::string& process(prov.processName());
0289         edm::InputTag tagNew(edm::InputTag(label, instance, process));
0290         if (tagOld.encode() != tagNew.encode()) {
0291           filterproduct.addCollectionTag(tagNew);
0292           tagOld = tagNew;
0293         }
0294       }
0295       for (unsigned int i = 0; i < originTag2_.size(); ++i) {
0296         filterproduct.addCollectionTag(originTag2_[i]);
0297       }
0298       tagOld = edm::InputTag();
0299       for (trigger::size_type i2 = 0; i2 != n2; ++i2) {
0300         const edm::ProductID pid(coll2[i2].id());
0301         const auto& prov = iEvent.getStableProvenance(pid);
0302         const std::string& label(prov.moduleLabel());
0303         const std::string& instance(prov.productInstanceName());
0304         const std::string& process(prov.processName());
0305         edm::InputTag tagNew(edm::InputTag(label, instance, process));
0306         if (tagOld.encode() != tagNew.encode()) {
0307           filterproduct.addCollectionTag(tagNew);
0308           tagOld = tagNew;
0309         }
0310       }
0311     }
0312 
0313     return true;
0314   } else
0315     return false;
0316 }
0317 
0318 template <typename T1, typename T2>
0319 bool HLTDoubletDZ<T1, T2>::computeDZ(edm::Event& iEvent, T1Ref& r1, T2Ref& r2) const {
0320   const reco::Candidate& candidate1(*r1);
0321   const reco::Candidate& candidate2(*r2);
0322   if (reco::deltaR(candidate1, candidate2) < minDR_)
0323     return false;
0324   if (std::abs(candidate1.vz() - candidate2.vz()) > maxDZ_)
0325     return false;
0326 
0327   return true;
0328 }
0329 
0330 template <>
0331 bool HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoChargedCandidate>::computeDZ(edm::Event& iEvent,
0332                                                                                   T1Ref& r1,
0333                                                                                   T2Ref& r2) const {
0334   edm::Handle<reco::ElectronCollection> electronHandle_;
0335   iEvent.getByToken(electronToken_, electronHandle_);
0336   if (!electronHandle_.isValid())
0337     edm::LogError("HLTDoubletDZ") << "HLTDoubletDZ: Electron Handle not valid.";
0338 
0339   if (reco::deltaR(*r1, *r2) < minDR_)
0340     return false;
0341   reco::Electron e1;
0342   for (auto const& eleIt : *electronHandle_) {
0343     if (eleIt.superCluster() == r1->superCluster())
0344       e1 = eleIt;
0345   }
0346 
0347   const reco::RecoChargedCandidate& candidate2(*r2);
0348   bool skipDZ = false;
0349   if (minPixHitsForDZ_ > 0 && (e1.gsfTrack()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_ ||
0350                                candidate2.track()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_))
0351     skipDZ = true;
0352   if (!skipDZ && std::abs(e1.vz() - candidate2.vz()) > maxDZ_)
0353     return false;
0354 
0355   return true;
0356 }
0357 
0358 template <>
0359 bool HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoEcalCandidate>::computeDZ(edm::Event& iEvent,
0360                                                                                   T1Ref& r1,
0361                                                                                   T2Ref& r2) const {
0362   edm::Handle<reco::ElectronCollection> electronHandle_;
0363   iEvent.getByToken(electronToken_, electronHandle_);
0364   if (!electronHandle_.isValid())
0365     edm::LogError("HLTDoubletDZ") << "HLTDoubletDZ: Electron Handle not valid.";
0366 
0367   if (reco::deltaR(*r1, *r2) < minDR_)
0368     return false;
0369   reco::Electron e2;
0370   for (auto const& eleIt : *electronHandle_) {
0371     if (eleIt.superCluster() == r2->superCluster())
0372       e2 = eleIt;
0373   }
0374 
0375   const reco::RecoChargedCandidate& candidate1(*r1);
0376   bool skipDZ = false;
0377   if (minPixHitsForDZ_ > 0 && (candidate1.track()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_ ||
0378                                e2.gsfTrack()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_))
0379     skipDZ = true;
0380   if (!skipDZ && std::abs(e2.vz() - candidate1.vz()) > maxDZ_)
0381     return false;
0382 
0383   return true;
0384 }
0385 
0386 template <>
0387 bool HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoEcalCandidate>::computeDZ(edm::Event& iEvent,
0388                                                                                T1Ref& r1,
0389                                                                                T2Ref& r2) const {
0390   edm::Handle<reco::ElectronCollection> electronHandle_;
0391   iEvent.getByToken(electronToken_, electronHandle_);
0392   if (!electronHandle_.isValid())
0393     edm::LogError("HLTDoubletDZ") << "HLTDoubletDZ: Electron Handle not valid.";
0394 
0395   if (reco::deltaR(*r1, *r2) < minDR_)
0396     return false;
0397   reco::Electron e1, e2;
0398   for (auto const& eleIt : *electronHandle_) {
0399     if (eleIt.superCluster() == r2->superCluster())
0400       e2 = eleIt;
0401     if (eleIt.superCluster() == r1->superCluster())
0402       e1 = eleIt;
0403   }
0404 
0405   bool skipDZ = false;
0406   if (minPixHitsForDZ_ > 0 && (e1.gsfTrack()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_ ||
0407                                e2.gsfTrack()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_))
0408     skipDZ = true;
0409   if (!skipDZ && std::abs(e2.vz() - e1.vz()) > maxDZ_)
0410     return false;
0411 
0412   return true;
0413 }
0414 
0415 template <>
0416 bool HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoChargedCandidate>::computeDZ(edm::Event& iEvent,
0417                                                                                      T1Ref& r1,
0418                                                                                      T2Ref& r2) const {
0419   const reco::RecoChargedCandidate& candidate1(*r1);
0420   const reco::RecoChargedCandidate& candidate2(*r2);
0421   if (reco::deltaR(candidate1, candidate2) < minDR_)
0422     return false;
0423   bool skipDZ = false;
0424   if (minPixHitsForDZ_ > 0 && (candidate1.track()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_ ||
0425                                candidate2.track()->hitPattern().numberOfValidPixelHits() < minPixHitsForDZ_))
0426     skipDZ = true;
0427   if (!skipDZ && std::abs(candidate1.vz() - candidate2.vz()) > maxDZ_)
0428     return false;
0429 
0430   return true;
0431 }
0432 
0433 template <>
0434 bool HLTDoubletDZ<l1t::TkMuon, l1t::TkMuon>::computeDZ(edm::Event& iEvent,
0435                                                        l1t::TkMuonRef& r1,
0436                                                        l1t::TkMuonRef& r2) const {
0437   const l1t::TkMuon& candidate1(*r1);
0438   const l1t::TkMuon& candidate2(*r2);
0439   if (reco::deltaR(candidate1, candidate2) < minDR_)
0440     return false;
0441 
0442   // We don't care about minPixHitsForDZ_ with the L1TkMuons,
0443   // especially because the pixel does not participate in the L1T
0444   if (std::abs(candidate1.trkzVtx() - candidate2.trkzVtx()) > maxDZ_)
0445     return false;
0446 
0447   return true;
0448 }
0449 
0450 template <>
0451 bool HLTDoubletDZ<l1t::HPSPFTau, l1t::HPSPFTau>::computeDZ(edm::Event& iEvent,
0452                                                            l1t::HPSPFTauRef& r1,
0453                                                            l1t::HPSPFTauRef& r2) const {
0454   const l1t::HPSPFTau& candidate1(*r1);
0455   const l1t::HPSPFTau& candidate2(*r2);
0456   if (reco::deltaR(candidate1, candidate2) < minDR_)
0457     return false;
0458 
0459   // We don't care about minPixHitsForDZ_ with the L1HPSPFTaus,
0460   // especially because the pixel does not participate in the L1T
0461   if (std::abs(candidate1.leadChargedPFCand()->pfTrack()->vertex().z() -
0462                candidate2.leadChargedPFCand()->pfTrack()->vertex().z()) > maxDZ_)
0463     return false;
0464 
0465   return true;
0466 }
0467 
0468 // ------------ method called to produce the data  ------------
0469 template <typename T1, typename T2>
0470 bool HLTDoubletDZ<T1, T2>::hltFilter(edm::Event& iEvent,
0471                                      const edm::EventSetup& iSetup,
0472                                      trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0473   // All HLT filters must create and fill an HLT filter object,
0474   // recording any reconstructed physics objects satisfying (or not)
0475   // this HLT filter, and place it in the Event.
0476   bool accept(false);
0477 
0478   std::vector<T1Ref> coll1;
0479   std::vector<T2Ref> coll2;
0480 
0481   if (getCollections(iEvent, coll1, coll2, filterproduct)) {
0482     int n(0);
0483     T1Ref r1;
0484     T2Ref r2;
0485 
0486     for (unsigned int i1 = 0; i1 != coll1.size(); i1++) {
0487       r1 = coll1[i1];
0488       //const reco::Candidate& candidate1(*r1);
0489       unsigned int I(0);
0490       if (same_) {
0491         I = i1 + 1;
0492       }
0493       for (unsigned int i2 = I; i2 != coll2.size(); i2++) {
0494         r2 = coll2[i2];
0495         if (checkSC_) {
0496           if (r1->superCluster().isNonnull() && r2->superCluster().isNonnull()) {
0497             if (r1->superCluster() == r2->superCluster())
0498               continue;
0499           }
0500         }
0501 
0502         if (!computeDZ(iEvent, r1, r2))
0503           continue;
0504 
0505         n++;
0506         filterproduct.addObject(triggerType1_, r1);
0507         filterproduct.addObject(triggerType2_, r2);
0508       }
0509     }
0510 
0511     accept = accept || (n >= min_N_);
0512   }
0513 
0514   return accept;
0515 }
0516 
0517 /// Special instantiation for L1TkMuon
0518 /// L1TkMuon are *not* RecoCandidates, therefore they don't implement superCluster()
0519 /// They are LeafCandidates instead
0520 template <>
0521 bool HLTDoubletDZ<l1t::TkMuon, l1t::TkMuon>::hltFilter(edm::Event& iEvent,
0522                                                        const edm::EventSetup& iSetup,
0523                                                        trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0524   // All HLT filters must create and fill an HLT filter object,
0525   // recording any reconstructed physics objects satisfying (or not)
0526   // this HLT filter, and place it in the Event.
0527   bool accept(false);
0528 
0529   std::vector<l1t::TkMuonRef> coll1;
0530   std::vector<l1t::TkMuonRef> coll2;
0531 
0532   if (getCollections(iEvent, coll1, coll2, filterproduct)) {
0533     int n(0);
0534     l1t::TkMuonRef r1;
0535     l1t::TkMuonRef r2;
0536 
0537     for (unsigned int i1 = 0; i1 != coll1.size(); i1++) {
0538       r1 = coll1[i1];
0539       unsigned int I(0);
0540       if (same_) {
0541         I = i1 + 1;
0542       }
0543       for (unsigned int i2 = I; i2 != coll2.size(); i2++) {
0544         r2 = coll2[i2];
0545 
0546         if (!computeDZ(iEvent, r1, r2))
0547           continue;
0548 
0549         n++;
0550         filterproduct.addObject(triggerType1_, r1);
0551         filterproduct.addObject(triggerType2_, r2);
0552       }
0553     }
0554 
0555     accept = accept || (n >= min_N_);
0556   }
0557 
0558   return accept;
0559 }
0560 
0561 /// Special instantiation for L1PFTau
0562 /// L1PFTau are *not* RecoCandidates, therefore they don't implement superCluster()
0563 /// They are LeafCandidates instead
0564 template <>
0565 bool HLTDoubletDZ<l1t::PFTau, l1t::PFTau>::hltFilter(edm::Event& iEvent,
0566                                                      const edm::EventSetup& iSetup,
0567                                                      trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0568   // All HLT filters must create and fill an HLT filter object,
0569   // recording any reconstructed physics objects satisfying (or not)
0570   // this HLT filter, and place it in the Event.
0571   bool accept(false);
0572 
0573   std::vector<l1t::PFTauRef> coll1;
0574   std::vector<l1t::PFTauRef> coll2;
0575 
0576   if (getCollections(iEvent, coll1, coll2, filterproduct)) {
0577     int n(0);
0578     l1t::PFTauRef r1;
0579     l1t::PFTauRef r2;
0580 
0581     for (unsigned int i1 = 0; i1 != coll1.size(); i1++) {
0582       r1 = coll1[i1];
0583       unsigned int I(0);
0584       if (same_) {
0585         I = i1 + 1;
0586       }
0587       for (unsigned int i2 = I; i2 != coll2.size(); i2++) {
0588         r2 = coll2[i2];
0589 
0590         if (!computeDZ(iEvent, r1, r2))
0591           continue;
0592 
0593         n++;
0594         filterproduct.addObject(triggerType1_, r1);
0595         filterproduct.addObject(triggerType2_, r2);
0596       }
0597     }
0598 
0599     accept = accept || (n >= min_N_);
0600   }
0601 
0602   return accept;
0603 }
0604 
0605 /// Special instantiation for L1HPSPFTau
0606 /// L1HPSPFTau are *not* RecoCandidates, therefore they don't implement superCluster()
0607 /// They are LeafCandidates instead
0608 template <>
0609 bool HLTDoubletDZ<l1t::HPSPFTau, l1t::HPSPFTau>::hltFilter(edm::Event& iEvent,
0610                                                            const edm::EventSetup& iSetup,
0611                                                            trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0612   // All HLT filters must create and fill an HLT filter object,
0613   // recording any reconstructed physics objects satisfying (or not)
0614   // this HLT filter, and place it in the Event.
0615   bool accept(false);
0616 
0617   std::vector<l1t::HPSPFTauRef> coll1;
0618   std::vector<l1t::HPSPFTauRef> coll2;
0619 
0620   if (getCollections(iEvent, coll1, coll2, filterproduct)) {
0621     int n(0);
0622     l1t::HPSPFTauRef r1;
0623     l1t::HPSPFTauRef r2;
0624 
0625     for (unsigned int i1 = 0; i1 != coll1.size(); i1++) {
0626       r1 = coll1[i1];
0627       unsigned int I(0);
0628       if (same_) {
0629         I = i1 + 1;
0630       }
0631       for (unsigned int i2 = I; i2 != coll2.size(); i2++) {
0632         r2 = coll2[i2];
0633 
0634         if (!computeDZ(iEvent, r1, r2))
0635           continue;
0636 
0637         n++;
0638         filterproduct.addObject(triggerType1_, r1);
0639         filterproduct.addObject(triggerType2_, r2);
0640       }
0641     }
0642 
0643     accept = accept || (n >= min_N_);
0644   }
0645 
0646   return accept;
0647 }
0648 
0649 typedef HLTDoubletDZ<reco::Electron, reco::Electron> HLT2ElectronElectronDZ;
0650 typedef HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoChargedCandidate> HLT2MuonMuonDZ;
0651 typedef HLTDoubletDZ<reco::Electron, reco::RecoChargedCandidate> HLT2ElectronMuonDZ;
0652 typedef HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoEcalCandidate> HLT2PhotonPhotonDZ;
0653 typedef HLTDoubletDZ<reco::RecoChargedCandidate, reco::RecoEcalCandidate> HLT2MuonPhotonDZ;
0654 typedef HLTDoubletDZ<reco::RecoEcalCandidate, reco::RecoChargedCandidate> HLT2PhotonMuonDZ;
0655 typedef HLTDoubletDZ<l1t::TkMuon, l1t::TkMuon> HLT2L1TkMuonL1TkMuonDZ;
0656 typedef HLTDoubletDZ<l1t::PFTau, l1t::PFTau> HLT2L1PFTauL1PFTauDZ;
0657 typedef HLTDoubletDZ<l1t::HPSPFTau, l1t::HPSPFTau> HLT2L1HPSPFTauL1HPSPFTauDZ;
0658 
0659 DEFINE_FWK_MODULE(HLT2ElectronElectronDZ);
0660 DEFINE_FWK_MODULE(HLT2MuonMuonDZ);
0661 DEFINE_FWK_MODULE(HLT2ElectronMuonDZ);
0662 DEFINE_FWK_MODULE(HLT2PhotonPhotonDZ);
0663 DEFINE_FWK_MODULE(HLT2PhotonMuonDZ);
0664 DEFINE_FWK_MODULE(HLT2MuonPhotonDZ);
0665 DEFINE_FWK_MODULE(HLT2L1TkMuonL1TkMuonDZ);
0666 DEFINE_FWK_MODULE(HLT2L1PFTauL1PFTauDZ);
0667 DEFINE_FWK_MODULE(HLT2L1HPSPFTauL1HPSPFTauDZ);