File indexing completed on 2023-03-17 11:09:26
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "HLTDoublet.h"
0012
0013 #include "DataFormats/Common/interface/Handle.h"
0014
0015 #include "DataFormats/Common/interface/Ref.h"
0016 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0017
0018 #include "DataFormats/Candidate/interface/Particle.h"
0019
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021
0022 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0023 #include <cmath>
0024
0025
0026
0027
0028 template <typename T1, typename T2>
0029 HLTDoublet<T1, T2>::HLTDoublet(const edm::ParameterSet& iConfig)
0030 : HLTFilter(iConfig),
0031 originTag1_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag1")),
0032 originTag2_(iConfig.template getParameter<std::vector<edm::InputTag>>("originTag2")),
0033 inputTag1_(iConfig.template getParameter<edm::InputTag>("inputTag1")),
0034 inputTag2_(iConfig.template getParameter<edm::InputTag>("inputTag2")),
0035 inputToken1_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag1_)),
0036 inputToken2_(consumes<trigger::TriggerFilterObjectWithRefs>(inputTag2_)),
0037 triggerType1_(iConfig.template getParameter<int>("triggerType1")),
0038 triggerType2_(iConfig.template getParameter<int>("triggerType2")),
0039 min_Dphi_(iConfig.template getParameter<double>("MinDphi")),
0040 max_Dphi_(iConfig.template getParameter<double>("MaxDphi")),
0041 min_Deta_(iConfig.template getParameter<double>("MinDeta")),
0042 max_Deta_(iConfig.template getParameter<double>("MaxDeta")),
0043 min_Minv_(iConfig.template getParameter<double>("MinMinv")),
0044 max_Minv_(iConfig.template getParameter<double>("MaxMinv")),
0045 min_DelR_(iConfig.template getParameter<double>("MinDelR")),
0046 max_DelR_(iConfig.template getParameter<double>("MaxDelR")),
0047 min_Pt_(iConfig.template getParameter<double>("MinPt")),
0048 max_Pt_(iConfig.template getParameter<double>("MaxPt")),
0049 min_N_(iConfig.template getParameter<int>("MinN")),
0050 same_(inputTag1_.encode() == inputTag2_.encode()),
0051 cutdphi_(min_Dphi_ <= max_Dphi_),
0052 cutdeta_(min_Deta_ <= max_Deta_),
0053 cutminv_(min_Minv_ <= max_Minv_),
0054 cutdelr_(min_DelR_ <= max_DelR_),
0055 cutpt_(min_Pt_ <= max_Pt_)
0056 {
0057 LogDebug("") << "InputTags and cuts : " << inputTag1_.encode() << " " << inputTag2_.encode() << triggerType1_ << " "
0058 << triggerType2_ << " Dphi [" << min_Dphi_ << " " << max_Dphi_ << "]"
0059 << " Deta [" << min_Deta_ << " " << max_Deta_ << "]"
0060 << " Minv [" << min_Minv_ << " " << max_Minv_ << "]"
0061 << " DelR [" << min_DelR_ << " " << max_DelR_ << "]"
0062 << " Pt [" << min_Pt_ << " " << max_Pt_ << "]"
0063 << " MinN =" << min_N_ << " same/dphi/deta/minv/delr/pt " << same_ << cutdphi_ << cutdeta_ << cutminv_
0064 << cutdelr_ << cutpt_;
0065 }
0066
0067 template <typename T1, typename T2>
0068 HLTDoublet<T1, T2>::~HLTDoublet() = default;
0069 template <typename T1, typename T2>
0070 void HLTDoublet<T1, T2>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0071 edm::ParameterSetDescription desc;
0072 makeHLTFilterDescription(desc);
0073 std::vector<edm::InputTag> originTag1(1, edm::InputTag("hltOriginal1"));
0074 std::vector<edm::InputTag> originTag2(1, edm::InputTag("hltOriginal2"));
0075 desc.add<std::vector<edm::InputTag>>("originTag1", originTag1);
0076 desc.add<std::vector<edm::InputTag>>("originTag2", originTag2);
0077 desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0078 desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered22"));
0079 desc.add<int>("triggerType1", 0);
0080 desc.add<int>("triggerType2", 0);
0081 desc.add<double>("MinDphi", +1.0);
0082 desc.add<double>("MaxDphi", -1.0);
0083 desc.add<double>("MinDeta", +1.0);
0084 desc.add<double>("MaxDeta", -1.0);
0085 desc.add<double>("MinMinv", +1.0);
0086 desc.add<double>("MaxMinv", -1.0);
0087 desc.add<double>("MinDelR", +1.0);
0088 desc.add<double>("MaxDelR", -1.0);
0089 desc.add<double>("MinPt", +1.0);
0090 desc.add<double>("MaxPt", -1.0);
0091 desc.add<int>("MinN", 1);
0092 descriptions.add(defaultModuleLabel<HLTDoublet<T1, T2>>(), desc);
0093 }
0094
0095
0096
0097
0098
0099
0100 template <typename T1, typename T2>
0101 bool HLTDoublet<T1, T2>::hltFilter(edm::Event& iEvent,
0102 const edm::EventSetup& iSetup,
0103 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0104 using namespace std;
0105 using namespace edm;
0106 using namespace reco;
0107 using namespace trigger;
0108
0109
0110
0111
0112
0113 bool accept(false);
0114
0115 LogVerbatim("HLTDoublet") << " XXX " << moduleLabel() << " 0 " << std::endl;
0116
0117 std::vector<T1Ref> coll1;
0118 std::vector<T2Ref> coll2;
0119
0120
0121 Handle<TriggerFilterObjectWithRefs> handle1, handle2;
0122 if (iEvent.getByToken(inputToken1_, handle1) && iEvent.getByToken(inputToken2_, handle2)) {
0123 handle1->getObjects(triggerType1_, coll1);
0124 handle2->getObjects(triggerType2_, coll2);
0125 const size_type n1(coll1.size());
0126 const size_type n2(coll2.size());
0127
0128 if (saveTags()) {
0129 InputTag tagOld;
0130 for (unsigned int i = 0; i < originTag1_.size(); ++i) {
0131 filterproduct.addCollectionTag(originTag1_[i]);
0132 LogVerbatim("HLTDoublet") << " XXX " << moduleLabel() << " 1a/" << i << " " << originTag1_[i].encode()
0133 << std::endl;
0134 }
0135 tagOld = InputTag();
0136 for (size_type i1 = 0; i1 != n1; ++i1) {
0137 const ProductID pid(coll1[i1].id());
0138 const auto& prov = iEvent.getStableProvenance(pid);
0139 const string& label(prov.moduleLabel());
0140 const string& instance(prov.productInstanceName());
0141 const string& process(prov.processName());
0142 InputTag tagNew(InputTag(label, instance, process));
0143 if (tagOld.encode() != tagNew.encode()) {
0144 filterproduct.addCollectionTag(tagNew);
0145 tagOld = tagNew;
0146 LogVerbatim("HLTDoublet") << " XXX " << moduleLabel() << " 1b " << tagNew.encode() << std::endl;
0147 }
0148 }
0149 for (unsigned int i = 0; i < originTag2_.size(); ++i) {
0150 filterproduct.addCollectionTag(originTag2_[i]);
0151 LogVerbatim("HLTDoublet") << " XXX " << moduleLabel() << " 2a/" << i << " " << originTag2_[i].encode()
0152 << std::endl;
0153 }
0154 tagOld = InputTag();
0155 for (size_type i2 = 0; i2 != n2; ++i2) {
0156 const ProductID pid(coll2[i2].id());
0157 const auto& prov = iEvent.getStableProvenance(pid);
0158 const string& label(prov.moduleLabel());
0159 const string& instance(prov.productInstanceName());
0160 const string& process(prov.processName());
0161 InputTag tagNew(InputTag(label, instance, process));
0162 if (tagOld.encode() != tagNew.encode()) {
0163 filterproduct.addCollectionTag(tagNew);
0164 tagOld = tagNew;
0165 LogVerbatim("HLTDoublet") << " XXX " << moduleLabel() << " 2b " << tagNew.encode() << std::endl;
0166 }
0167 }
0168 }
0169
0170 int n(0);
0171 T1Ref r1;
0172 T2Ref r2;
0173 Particle::LorentzVector p1, p2, p;
0174 for (unsigned int i1 = 0; i1 != n1; i1++) {
0175 r1 = coll1[i1];
0176 p1 = r1->p4();
0177 unsigned int I(0);
0178 if (same_) {
0179 I = i1 + 1;
0180 }
0181 for (unsigned int i2 = I; i2 != n2; i2++) {
0182 r2 = coll2[i2];
0183 p2 = r2->p4();
0184
0185 double Dphi(std::abs(p1.phi() - p2.phi()));
0186 if (Dphi > M_PI)
0187 Dphi = 2.0 * M_PI - Dphi;
0188
0189 double Deta(std::abs(p1.eta() - p2.eta()));
0190
0191 p = p1 + p2;
0192 double Minv(std::abs(p.mass()));
0193 double Pt(p.pt());
0194
0195 double DelR(sqrt(Dphi * Dphi + Deta * Deta));
0196
0197 if (((!cutdphi_) || ((min_Dphi_ <= Dphi) && (Dphi <= max_Dphi_))) &&
0198 ((!cutdeta_) || ((min_Deta_ <= Deta) && (Deta <= max_Deta_))) &&
0199 ((!cutminv_) || ((min_Minv_ <= Minv) && (Minv <= max_Minv_))) &&
0200 ((!cutdelr_) || ((min_DelR_ <= DelR) && (DelR <= max_DelR_))) &&
0201 ((!cutpt_) || ((min_Pt_ <= Pt) && (Pt <= max_Pt_)))) {
0202 n++;
0203 filterproduct.addObject(triggerType1_, r1);
0204 filterproduct.addObject(triggerType2_, r2);
0205 }
0206 }
0207 }
0208
0209 accept = (n >= min_N_);
0210 }
0211
0212 return accept;
0213 }