HLTfilters/plugins/HLTDoubletSinglet.h

0001 #ifndef HLTrigger_HLTfilters_plugins_HLTDoubletSinglet_h
0002 #define HLTrigger_HLTfilters_plugins_HLTDoubletSinglet_h
0003
0004 /** \class HLTDoubletSinglet
0005  *
0006  *
0007  *  This class is an HLTFilter (-> EDFilter) implementing a basic HLT
0008  *  trigger for triplets of objects, evaluating all triplets with the first
0009  *  object from collection 1, the second object from collection 2,
0010  *  and the third object from collection 3,
0011  *  cutting on variables relating to their 4-momentum representations.
0012  *  The filter itself compares only objects from collection 3
0013  *  with objects from collections 1 and 2.
0014  *  The object collections are assumed to be outputs of HLTSinglet
0015  *  single-object-type filters so that the access is thorugh
0016  *  RefToBases and polymorphic.
0017  *
0018  *
0019  *  \author Jaime Leon Holgado
0020  *
0021  */
0022
0023 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0024 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0025
0026 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0027 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0028
0029 #include "DataFormats/Common/interface/Ref.h"
0030 #include "DataFormats/Common/interface/Handle.h"
0031 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0032 #include "DataFormats/Candidate/interface/Candidate.h"
0033 #include "DataFormats/Math/interface/deltaPhi.h"
0034
0035 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0036
0037 #include <string>
0038 #include <vector>
0039 #include <cmath>
0040
0041 //
0042 // class declaration
0043 //
0044
0045 template <typename T1, typename T2, typename T3>
0046 class HLTDoubletSinglet : public HLTFilter {
0047 public:
0048   explicit HLTDoubletSinglet(const edm::ParameterSet&);
0049   ~HLTDoubletSinglet() override;
0050   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0051   bool hltFilter(edm::Event&,
0052                  const edm::EventSetup&,
0053                  trigger::TriggerFilterObjectWithRefs& filterproduct) const override;
0054
0055 private:
0056   // configuration
0057   const std::vector<edm::InputTag> originTag1_;  // input tag identifying originals 1st product
0058   const std::vector<edm::InputTag> originTag2_;  // input tag identifying originals 2nd product
0059   const std::vector<edm::InputTag> originTag3_;  // input tag identifying originals 3rd product
0060   const edm::InputTag inputTag1_;                // input tag identifying filtered 1st product
0061   const edm::InputTag inputTag2_;                // input tag identifying filtered 2nd product
0062   const edm::InputTag inputTag3_;                // input tag identifying filtered 3rd product
0063   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> inputToken1_;
0064   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> inputToken2_;
0065   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> inputToken3_;
0066   const int triggerType1_;
0067   const int triggerType2_;
0068   const int triggerType3_;
0069   const double min_Dphi_, max_Dphi_;  // DeltaPhi (1,3) and (2,3) window
0070   const double min_Deta_, max_Deta_;  // DeltaEta (1,3) and (2,3) window
0071   const double min_Minv_, max_Minv_;  // Minv(1,2) and Minv(2,3) window
0072   const double min_DelR_, max_DelR_;  // DeltaR (1,3) and (2,3) window
0073   const double min_Pt_, max_Pt_;      // Pt(1,3) and (2,3) window
0074   const int min_N_;                   // number of triplets passing cuts required
0075
0076   // calculated from configuration in c'tor
0077   const bool same12_, same13_, same23_;                       // 1st and 2nd product are one and the same
0078   const double min_DelR2_, max_DelR2_;                        // DeltaR (1,3) and (2,3) window
0079   const bool cutdphi_, cutdeta_, cutminv_, cutdelr_, cutpt_;  // cuts are on=true or off=false
0080
0081   // typedefs
0082   typedef std::vector<T1> T1Collection;
0083   typedef edm::Ref<T1Collection> T1Ref;
0084   typedef std::vector<T2> T2Collection;
0085   typedef edm::Ref<T2Collection> T2Ref;
0086   typedef std::vector<T3> T3Collection;
0087   typedef edm::Ref<T3Collection> T3Ref;
0088 };
0089
0090 //
0091 // class implementation
0092 //
0093
0094 //
0095 // constructors and destructor
0096 //
0097 template <typename T1, typename T2, typename T3>
0098 HLTDoubletSinglet<T1, T2, T3>::HLTDoubletSinglet(const edm::ParameterSet& iConfig)
0099     : HLTFilter(iConfig),
0100       originTag1_(iConfig.getParameter<std::vector<edm::InputTag>>("originTag1")),
0101       originTag2_(iConfig.getParameter<std::vector<edm::InputTag>>("originTag2")),
0102       originTag3_(iConfig.getParameter<std::vector<edm::InputTag>>("originTag3")),
0103       inputTag1_(iConfig.getParameter<edm::InputTag>("inputTag1")),
0104       inputTag2_(iConfig.getParameter<edm::InputTag>("inputTag2")),
0105       inputTag3_(iConfig.getParameter<edm::InputTag>("inputTag3")),
0106       inputToken1_(consumes(inputTag1_)),
0107       inputToken2_(consumes(inputTag2_)),
0108       inputToken3_(consumes(inputTag3_)),
0109       triggerType1_(iConfig.getParameter<int>("triggerType1")),
0110       triggerType2_(iConfig.getParameter<int>("triggerType2")),
0111       triggerType3_(iConfig.getParameter<int>("triggerType3")),
0112       min_Dphi_(iConfig.getParameter<double>("MinDphi")),
0113       max_Dphi_(iConfig.getParameter<double>("MaxDphi")),
0114       min_Deta_(iConfig.getParameter<double>("MinDeta")),
0115       max_Deta_(iConfig.getParameter<double>("MaxDeta")),
0116       min_Minv_(iConfig.getParameter<double>("MinMinv")),
0117       max_Minv_(iConfig.getParameter<double>("MaxMinv")),
0118       min_DelR_(iConfig.getParameter<double>("MinDelR")),
0119       max_DelR_(iConfig.getParameter<double>("MaxDelR")),
0120       min_Pt_(iConfig.getParameter<double>("MinPt")),
0121       max_Pt_(iConfig.getParameter<double>("MaxPt")),
0122       min_N_(iConfig.getParameter<int>("MinN")),
0123       same12_(inputTag1_.encode() == inputTag2_.encode()),    // same collections to be compared?
0124       same13_(inputTag1_.encode() == inputTag3_.encode()),    // same collections to be compared?
0125       same23_(inputTag2_.encode() == inputTag3_.encode()),    // same collections to be compared?
0126       min_DelR2_(min_DelR_ < 0 ? 0 : min_DelR_ * min_DelR_),  // avoid computing sqrt(R2)
0127       max_DelR2_(max_DelR_ < 0 ? 0 : max_DelR_ * max_DelR_),  // avoid computing sqrt(R2)
0128       cutdphi_(min_Dphi_ <= max_Dphi_),                       // cut active?
0129       cutdeta_(min_Deta_ <= max_Deta_),                       // cut active?
0130       cutminv_(min_Minv_ <= max_Minv_),                       // cut active?
0131       cutdelr_(min_DelR_ <= max_DelR_),                       // cut active?
0132       cutpt_(min_Pt_ <= max_Pt_)                              // cut active?
0133 {
0134   LogDebug("") << "InputTags and cuts : " << inputTag1_.encode() << " " << inputTag2_.encode() << " "
0135                << inputTag3_.encode() << triggerType1_ << " " << triggerType2_ << " " << triggerType3_ << " Dphi ["
0136                << min_Dphi_ << " " << max_Dphi_ << "]"
0137                << " Deta [" << min_Deta_ << " " << max_Deta_ << "]"
0138                << " Minv [" << min_Minv_ << " " << max_Minv_ << "]"
0139                << " DelR [" << min_DelR_ << " " << max_DelR_ << "]"
0140                << " Pt   [" << min_Pt_ << " " << max_Pt_ << "]"
0141                << " MinN =" << min_N_ << " same12/same13/same23/dphi/deta/minv/delr/pt " << same12_ << same13_
0142                << same23_ << cutdphi_ << cutdeta_ << cutminv_ << cutdelr_ << cutpt_;
0143   if (cutdelr_ && max_DelR_ <= 0)
0144     edm::LogWarning("HLTDoubletSinglet")
0145         << " moduleLabel: " << moduleLabel()
0146         << "Warning: The deltaR requirement is active, but its range is invalid: DelR [" << min_DelR_ << " "
0147         << max_DelR_ << "]";
0148 }
0149
0150 template <typename T1, typename T2, typename T3>
0151 HLTDoubletSinglet<T1, T2, T3>::~HLTDoubletSinglet() = default;
0152
0153 template <typename T1, typename T2, typename T3>
0154 void HLTDoubletSinglet<T1, T2, T3>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0155   edm::ParameterSetDescription desc;
0156   makeHLTFilterDescription(desc);
0157   desc.add<std::vector<edm::InputTag>>("originTag1", {edm::InputTag("hltOriginal1")});
0158   desc.add<std::vector<edm::InputTag>>("originTag2", {edm::InputTag("hltOriginal2")});
0159   desc.add<std::vector<edm::InputTag>>("originTag3", {edm::InputTag("hltOriginal3")});
0160   desc.add<edm::InputTag>("inputTag1", edm::InputTag("hltFiltered1"));
0161   desc.add<edm::InputTag>("inputTag2", edm::InputTag("hltFiltered2"));
0162   desc.add<edm::InputTag>("inputTag3", edm::InputTag("hltFiltered3"));
0163   desc.add<int>("triggerType1", 0);
0164   desc.add<int>("triggerType2", 0);
0165   desc.add<int>("triggerType3", 0);
0166   desc.add<double>("MinDphi", +1.0);
0167   desc.add<double>("MaxDphi", -1.0);
0168   desc.add<double>("MinDeta", +1.0);
0169   desc.add<double>("MaxDeta", -1.0);
0170   desc.add<double>("MinMinv", +1.0);
0171   desc.add<double>("MaxMinv", -1.0);
0172   desc.add<double>("MinDelR", +1.0);
0173   desc.add<double>("MaxDelR", -1.0);
0174   desc.add<double>("MinPt", +1.0);
0175   desc.add<double>("MaxPt", -1.0);
0176   desc.add<int>("MinN", 1);
0177   descriptions.add(defaultModuleLabel<HLTDoubletSinglet<T1, T2, T3>>(), desc);
0178 }
0179
0180 //
0181 // member functions
0182 //
0183
0184 // ------------ method called to produce the data  ------------
0185 template <typename T1, typename T2, typename T3>
0186 bool HLTDoubletSinglet<T1, T2, T3>::hltFilter(edm::Event& iEvent,
0187                                               const edm::EventSetup& iSetup,
0188                                               trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0189   using namespace std;
0190   using namespace edm;
0191   using namespace reco;
0192   using namespace trigger;
0193
0194   // All HLT filters must create and fill an HLT filter object,
0195   // recording any reconstructed physics objects satisfying (or not)
0196   // this HLT filter, and place it in the Event.
0197
0198   int n(0);
0199
0200   LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 0 ";
0201
0202   // get hold of pre-filtered object collections
0203   std::vector<T1Ref> coll1;
0204   auto const& objsWithRefs1 = iEvent.get(inputToken1_);
0205   objsWithRefs1.getObjects(triggerType1_, coll1);
0206   std::vector<T2Ref> coll2;
0207   auto const& objsWithRefs2 = iEvent.get(inputToken2_);
0208   objsWithRefs2.getObjects(triggerType2_, coll2);
0209   std::vector<T3Ref> coll3;
0210   auto const& objsWithRefs3 = iEvent.get(inputToken3_);
0211   objsWithRefs3.getObjects(triggerType3_, coll3);
0212
0213   const size_type n1(coll1.size());
0214   const size_type n2(coll2.size());
0215   const size_type n3(coll3.size());
0216
0217   if (saveTags()) {
0218     InputTag tagOld;
0219     for (size_t i = 0; i < originTag1_.size(); ++i) {
0220       filterproduct.addCollectionTag(originTag1_[i]);
0221       LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 1a/" << i << " "
0222                                        << originTag1_[i].encode();
0223     }
0224     tagOld = InputTag();
0225     for (size_type i1 = 0; i1 != n1; ++i1) {
0226       const ProductID pid(coll1[i1].id());
0227       const auto& prov = iEvent.getStableProvenance(pid);
0228       const string& label(prov.moduleLabel());
0229       const string& instance(prov.productInstanceName());
0230       const string& process(prov.processName());
0231       InputTag tagNew(InputTag(label, instance, process));
0232       if (tagOld.encode() != tagNew.encode()) {
0233         filterproduct.addCollectionTag(tagNew);
0234         tagOld = tagNew;
0235         LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 1b " << tagNew.encode();
0236       }
0237     }
0238     for (size_t i = 0; i < originTag2_.size(); ++i) {
0239       filterproduct.addCollectionTag(originTag2_[i]);
0240       LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 2a/" << i << " "
0241                                        << originTag2_[i].encode();
0242     }
0243     tagOld = InputTag();
0244     for (size_type i2 = 0; i2 != n2; ++i2) {
0245       const ProductID pid(coll2[i2].id());
0246       const auto& prov = iEvent.getStableProvenance(pid);
0247       const string& label(prov.moduleLabel());
0248       const string& instance(prov.productInstanceName());
0249       const string& process(prov.processName());
0250       InputTag tagNew(InputTag(label, instance, process));
0251       if (tagOld.encode() != tagNew.encode()) {
0252         filterproduct.addCollectionTag(tagNew);
0253         tagOld = tagNew;
0254         LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 2b " << tagNew.encode();
0255       }
0256     }
0257     for (size_t i = 0; i < originTag3_.size(); ++i) {
0258       filterproduct.addCollectionTag(originTag3_[i]);
0259       LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 3a/" << i << " "
0260                                        << originTag3_[i].encode();
0261     }
0262     tagOld = InputTag();
0263     for (size_type i3 = 0; i3 != n3; ++i3) {
0264       const ProductID pid(coll3[i3].id());
0265       const auto& prov = iEvent.getStableProvenance(pid);
0266       const string& label(prov.moduleLabel());
0267       const string& instance(prov.productInstanceName());
0268       const string& process(prov.processName());
0269       InputTag tagNew(InputTag(label, instance, process));
0270       if (tagOld.encode() != tagNew.encode()) {
0271         filterproduct.addCollectionTag(tagNew);
0272         tagOld = tagNew;
0273         LogVerbatim("HLTDoubletSinglet") << " moduleLabel: " << moduleLabel() << " 3b " << tagNew.encode();
0274       }
0275     }
0276
0277     T1Ref r1;
0278     T2Ref r2;
0279     T3Ref r3;
0280     Candidate::LorentzVector p1, p2, p3, p13, p23;
0281     for (size_t i1 = 0; i1 != n1; i1++) {
0282       r1 = coll1[i1];
0283       p1 = r1->p4();
0284       auto const i2_min = (same12_ ? i1 + 1 : 0);
0285       for (size_t i2 = i2_min; i2 != n2; i2++) {
0286         r2 = coll2[i2];
0287         p2 = r2->p4();
0288
0289         auto const i3_min = (same23_ ? i2_min + 1 : (same13_ ? i1 + 1 : 0));
0290         for (size_t i3 = i3_min; i3 != n3; i3++) {
0291           r3 = coll3[i3];
0292           p3 = r3->p4();
0293
0294           //deltaPhi
0295           auto const dPhi13(std::abs(deltaPhi(p1.phi(), p3.phi())));
0296           if (cutdphi_ && (min_Dphi_ > dPhi13 || dPhi13 > max_Dphi_))
0297             continue;
0298           auto const dPhi23(std::abs(deltaPhi(p2.phi(), p3.phi())));
0299           if (cutdphi_ && (min_Dphi_ > dPhi23 || dPhi23 > max_Dphi_))
0300             continue;
0301
0302           //deltaEta
0303           auto const dEta13(std::abs(p1.eta() - p3.eta()));
0304           if (cutdeta_ && (min_Deta_ > dEta13 || dEta13 > max_Deta_))
0305             continue;
0306           auto const dEta23(std::abs(p2.eta() - p3.eta()));
0307           if (cutdeta_ && (min_Deta_ > dEta23 || dEta23 > max_Deta_))
0308             continue;
0309
0310           //deltaR
0311           auto const delR2_13(dPhi13 * dPhi13 + dEta13 * dEta13);
0312           if (cutdelr_ && (min_DelR2_ > delR2_13 || delR2_13 > max_DelR2_))
0313             continue;
0314           auto const delR2_23(dPhi23 * dPhi23 + dEta23 * dEta23);
0315           if (cutdelr_ && (min_DelR2_ > delR2_23 || delR2_23 > max_DelR2_))
0316             continue;
0317
0318           //Pt and Minv
0319           p13 = p1 + p3;
0320           auto const mInv13(std::abs(p13.mass()));
0321           if (cutminv_ && (min_Minv_ > mInv13 || mInv13 > max_Minv_))
0322             continue;
0323           auto const pt13(p13.pt());
0324           if (cutpt_ && (min_Pt_ > pt13 || pt13 > max_Pt_))
0325             continue;
0326
0327           p23 = p2 + p3;
0328           auto const mInv23(std::abs(p23.mass()));
0329           if (cutminv_ && (min_Minv_ > mInv23 || mInv23 > max_Minv_))
0330             continue;
0331           auto const pt23(p23.pt());
0332           if (cutpt_ && (min_Pt_ > pt23 || pt23 > max_Pt_))
0333             continue;
0334
0335           n++;
0336           filterproduct.addObject(triggerType1_, r1);
0337           filterproduct.addObject(triggerType2_, r2);
0338           filterproduct.addObject(triggerType3_, r3);
0339         }
0340       }
0341     }
0342   }
0343
0344   // filter decision
0345   return (n >= min_N_);
0346 }
0347
0348 #endif  //HLTrigger_HLTfilters_plugins_HLTDoubletSinglet_h