Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:18:34 

0001 /** \class HLTJetSortedVBFFilter
0002  *
0003  * See header file for documentation
0004  *
0005  *
0006  *  \author Jacopo Bernardini
0007  *
0008  */
0009 #include <cmath>
0010 #include <string>
0011 #include <vector>
0012 
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016 #include "DataFormats/Common/interface/Handle.h"
0017 #include "DataFormats/Math/interface/deltaPhi.h"
0018 #include "DataFormats/Math/interface/deltaR.h"
0019 #include "DataFormats/Common/interface/RefToBase.h"
0020 #include "HLTrigger/JetMET/interface/HLTJetSortedVBFFilter.h"
0021 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0022 
0023 using namespace std;
0024 //
0025 // constructors and destructor//
0026 //
0027 template <typename T>
0028 HLTJetSortedVBFFilter<T>::HLTJetSortedVBFFilter(const edm::ParameterSet& iConfig)
0029     : HLTFilter(iConfig),
0030       inputJets_(iConfig.getParameter<edm::InputTag>("inputJets")),
0031       inputJetTags_(iConfig.getParameter<edm::InputTag>("inputJetTags")),
0032       mqq_(iConfig.getParameter<double>("Mqq")),
0033       detaqq_(iConfig.getParameter<double>("Detaqq")),
0034       detabb_(iConfig.getParameter<double>("Detabb")),
0035       dphibb_(iConfig.getParameter<double>("Dphibb")),
0036       ptsqq_(iConfig.getParameter<double>("Ptsumqq")),
0037       ptsbb_(iConfig.getParameter<double>("Ptsumbb")),
0038       seta_(iConfig.getParameter<double>("Etaq1Etaq2")),
0039       njets_(iConfig.getParameter<int>("njets")),
0040       value_(iConfig.getParameter<std::string>("value")),
0041       triggerType_(iConfig.getParameter<int>("triggerType")) {
0042   m_theJetsToken = consumes<std::vector<T>>(inputJets_);
0043   m_theJetTagsToken = consumes<reco::JetTagCollection>(inputJetTags_);
0044   if (njets_ < 4) {
0045     edm::LogWarning("LowNJets") << "njets=" << njets_ << " it must be >=4. Forced njets=4.";
0046     njets_ = 4;
0047   }
0048 }
0049 
0050 template <typename T>
0051 void HLTJetSortedVBFFilter<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0052   edm::ParameterSetDescription desc;
0053   makeHLTFilterDescription(desc);
0054   desc.add<edm::InputTag>("inputJets", edm::InputTag("hltJetCollection"));
0055   desc.add<edm::InputTag>("inputJetTags", edm::InputTag(""));
0056   desc.add<double>("Mqq", 200);
0057   desc.add<double>("Detaqq", 2.5);
0058   desc.add<double>("Detabb", 10.);
0059   desc.add<double>("Dphibb", 10.);
0060   desc.add<double>("Ptsumqq", 0.);
0061   desc.add<double>("Ptsumbb", 0.);
0062   desc.add<double>("Etaq1Etaq2", 40.);
0063   desc.add<std::string>("value", "second");
0064   desc.add<int>("triggerType", trigger::TriggerJet);
0065   desc.add<int>("njets", 4);
0066   descriptions.add(defaultModuleLabel<HLTJetSortedVBFFilter<T>>(), desc);
0067 }
0068 
0069 template <typename T>
0070 float HLTJetSortedVBFFilter<T>::findCSV(const typename std::vector<T>::const_iterator& jet,
0071                                         const reco::JetTagCollection& jetTags) {
0072   float minDr2 = 0.01f;
0073   float tmpCSV = -20;
0074   for (auto jetb = jetTags.begin(); (jetb != jetTags.end()); ++jetb) {
0075     float tmpDr2 = reco::deltaR2(*jet, *(jetb->first));
0076     if (tmpDr2 < minDr2) {
0077       minDr2 = tmpDr2;
0078       tmpCSV = jetb->second;
0079     }
0080   }
0081   return tmpCSV;
0082 }
0083 //
0084 // member functions
0085 //
0086 
0087 // ------------ method called to produce the data  ------------
0088 template <typename T>
0089 bool HLTJetSortedVBFFilter<T>::hltFilter(edm::Event& event,
0090                                          const edm::EventSetup& setup,
0091                                          trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0092   using namespace std;
0093   using namespace edm;
0094   using namespace reco;
0095   using namespace trigger;
0096 
0097   typedef vector<T> TCollection;
0098   typedef Ref<TCollection> TRef;
0099 
0100   bool accept(false);
0101 
0102   Handle<TCollection> jets;
0103   event.getByToken(m_theJetsToken, jets);
0104   Handle<JetTagCollection> jetTags;
0105 
0106   if (saveTags())
0107     filterproduct.addCollectionTag(inputJets_);
0108   if (jets->size() < 4)
0109     return false;
0110 
0111   const unsigned int nMax(njets_ < jets->size() ? njets_ : jets->size());
0112   vector<Jpair> sorted(nMax);
0113   vector<TRef> jetRefs(nMax);
0114   unsigned int nJetRefs = 0;
0115 
0116   unsigned int nJet = 0;
0117   double value(0.0);
0118 
0119   Particle::LorentzVector b1, b2, q1, q2;
0120   if (inputJetTags_.encode().empty()) {
0121     for (typename TCollection::const_iterator jet = jets->begin(); (jet != jets->end() && nJet < nMax); ++jet) {
0122       if (value_ == "Pt") {
0123         value = jet->pt();
0124       } else if (value_ == "Eta") {
0125         value = jet->eta();
0126       } else if (value_ == "Phi") {
0127         value = jet->phi();
0128       } else {
0129         value = 0.0;
0130       }
0131       sorted[nJet] = make_pair(value, nJet);
0132       ++nJet;
0133     }
0134     sort(sorted.begin(), sorted.end(), comparator);
0135     for (unsigned int i = 0; i < nMax; ++i) {
0136       jetRefs[i] = TRef(jets, sorted[i].second);
0137     }
0138     nJetRefs = nMax;
0139     q1 = jetRefs[3]->p4();
0140     b1 = jetRefs[2]->p4();
0141     b2 = jetRefs[1]->p4();
0142     q2 = jetRefs[0]->p4();
0143   } else if (value_ == "1BTagAndEta") {
0144     event.getByToken(m_theJetTagsToken, jetTags);
0145     vector<Jpair> sorted;
0146     unsigned int b1_idx = -1;
0147     float csv_max = -999;
0148     for (typename TCollection::const_iterator jet = jets->begin(); (jet != jets->end() && nJet < nMax);
0149          ++jet) {  //fill "sorted" and get the most b-tagged jet with higher CSV (b1)
0150       value = findCSV(jet, *jetTags);
0151       if (value > csv_max) {
0152         csv_max = value;
0153         b1_idx = nJet;
0154       }
0155       sorted.push_back(make_pair(jet->eta(), nJet));
0156       nJet++;
0157       //        cout << "jetPt=" << jet->pt() << "\tjetEta=" << jet->eta() << "\tjetCSV=" << value << endl;
0158     }
0159     if (b1_idx >= sorted.size())
0160       edm::LogError("OutOfRange") << "b1 index out of range.";
0161     sorted.erase(sorted.begin() + b1_idx);           //remove the most b-tagged jet from "sorted"
0162     sort(sorted.begin(), sorted.end(), comparator);  //sort "sorted" by eta
0163 
0164     unsigned int q1_idx = sorted.front().second;  //take the backward jet (q1)
0165     unsigned int q2_idx = sorted.back().second;   //take the forward jet (q2)
0166 
0167     unsigned int i = 0;
0168     while ((i == q1_idx) || (i == q2_idx) || (i == b1_idx))
0169       i++;  //take jet with highest pT but q1,q2,b1 (q2)
0170     unsigned int b2_idx = i;
0171 
0172     if (q1_idx < jets->size())
0173       q1 = jets->at(q1_idx).p4();
0174     else
0175       edm::LogWarning("Something wrong with q1");
0176     if (q2_idx < jets->size())
0177       q2 = jets->at(q2_idx).p4();
0178     else
0179       edm::LogWarning("Something wrong with q2");
0180     if (b1_idx < jets->size())
0181       b1 = jets->at(b1_idx).p4();
0182     else
0183       edm::LogWarning("Something wrong with b1");
0184     if (b2_idx < jets->size())
0185       b2 = jets->at(b2_idx).p4();
0186     else
0187       edm::LogWarning("Something wrong with b2");
0188 
0189     jetRefs[0] = TRef(jets, b1_idx);
0190     jetRefs[1] = TRef(jets, b2_idx);
0191     jetRefs[2] = TRef(jets, q1_idx);
0192     jetRefs[3] = TRef(jets, q2_idx);
0193     nJetRefs = 4;
0194 
0195     //      cout<<"\tPathB: b1="<<b1.pt()<<" b2="<<b2.pt()<<" q1="<<q1.pt()<<" q2="<<q2.pt()<<endl;
0196   } else if (value_ == "2BTagAndPt") {
0197     event.getByToken(m_theJetTagsToken, jetTags);
0198     vector<Jpair> sorted;
0199 
0200     unsigned int b1_idx = -1;
0201     unsigned int b2_idx = -1;
0202     float csv1 = -999;
0203     float csv2 = -999;
0204     for (typename TCollection::const_iterator jet = jets->begin(); (jet != jets->end() && nJet < nMax);
0205          ++jet) {  //fill "sorted" and get the two most b-tagged jets (b1,b2)
0206       value = findCSV(jet, *jetTags);
0207       if (value > csv1) {
0208         csv2 = csv1;
0209         b2_idx = b1_idx;
0210         csv1 = value;
0211         b1_idx = nJet;
0212       } else if (value > csv2) {
0213         csv2 = value;
0214         b2_idx = nJet;
0215       }
0216       sorted.push_back(make_pair(jet->eta(), nJet));
0217       nJet++;
0218       //        cout << "jetPt=" << jet->pt() << "\tjetEta=" << jet->eta() << "\tjetCSV=" << value << endl;
0219     }
0220     sorted.erase(sorted.begin() + b1_idx);  //remove b1 and b2 from sorted
0221     sorted.erase(sorted.begin() + (b1_idx > b2_idx ? b2_idx : b2_idx - 1));
0222 
0223     unsigned int q1_idx = sorted.at(0).second;  //get q1 and q2 as the jets with highest pT, but b1 and b2.
0224     unsigned int q2_idx = sorted.at(1).second;
0225 
0226     if (q1_idx < jets->size())
0227       q1 = jets->at(q1_idx).p4();
0228     else
0229       edm::LogWarning("Something wrong with q1");
0230     if (q2_idx < jets->size())
0231       q2 = jets->at(q2_idx).p4();
0232     else
0233       edm::LogWarning("Something wrong with q2");
0234     if (b1_idx < jets->size())
0235       b1 = jets->at(b1_idx).p4();
0236     else
0237       edm::LogWarning("Something wrong with b1");
0238     if (b2_idx < jets->size())
0239       b2 = jets->at(b2_idx).p4();
0240     else
0241       edm::LogWarning("Something wrong with b2");
0242 
0243     jetRefs[0] = TRef(jets, b1_idx);
0244     jetRefs[1] = TRef(jets, b2_idx);
0245     jetRefs[2] = TRef(jets, q1_idx);
0246     jetRefs[3] = TRef(jets, q2_idx);
0247     nJetRefs = 4;
0248 
0249     //      cout<<"\tPathA: b1="<<b1.pt()<<" b2="<<b2.pt()<<" q1="<<q1.pt()<<" q2="<<q2.pt()<<endl;
0250   } else {
0251     event.getByToken(m_theJetTagsToken, jetTags);
0252     for (typename TCollection::const_iterator jet = jets->begin(); (jet != jets->end() && nJet < nMax); ++jet) {
0253       if (value_ == "second") {
0254         value = findCSV(jet, *jetTags);
0255       } else {
0256         value = 0.0;
0257       }
0258       sorted[nJet] = make_pair(value, nJet);
0259       ++nJet;
0260     }
0261     sort(sorted.begin(), sorted.end(), comparator);
0262     for (unsigned int i = 0; i < nMax; ++i) {
0263       jetRefs[i] = TRef(jets, sorted[i].second);
0264     }
0265     nJetRefs = nMax;
0266     b1 = jetRefs[3]->p4();
0267     b2 = jetRefs[2]->p4();
0268     q1 = jetRefs[1]->p4();
0269     q2 = jetRefs[0]->p4();
0270   }
0271 
0272   double mqq_bs = (q1 + q2).M();
0273   double deltaetaqq = std::abs(q1.Eta() - q2.Eta());
0274   double deltaetabb = std::abs(b1.Eta() - b2.Eta());
0275   double deltaphibb = std::abs(reco::deltaPhi(b1.Phi(), b2.Phi()));
0276   double ptsqq_bs = (q1 + q2).Pt();
0277   double ptsbb_bs = (b1 + b2).Pt();
0278   double signeta = q1.Eta() * q2.Eta();
0279 
0280   if ((mqq_bs > mqq_) && (deltaetaqq > detaqq_) && (deltaetabb < detabb_) && (deltaphibb < dphibb_) &&
0281       (ptsqq_bs > ptsqq_) && (ptsbb_bs > ptsbb_) && (signeta < seta_)) {
0282     accept = true;
0283     for (unsigned int i = 0; i < nJetRefs; ++i) {
0284       filterproduct.addObject(triggerType_, jetRefs[i]);
0285     }
0286   }
0287 
0288   return accept;
0289 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!