HLTfilters/plugins/HLTDoublet.cc

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