File indexing completed on 2024-04-06 12:18:30
0001 #ifndef HLTCAWZTagFilter_h
0002 #define HLTCAWZTagFilter_h
0003
0004
0005 #include <memory>
0006 #include <vector>
0007 #include <sstream>
0008
0009
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "DataFormats/BTauReco/interface/CATopJetTagInfo.h"
0013 #include "FWCore/Utilities/interface/InputTag.h"
0014 #include "DataFormats/JetReco/interface/BasicJet.h"
0015 #include "DataFormats/JetReco/interface/CaloJet.h"
0016 #include "DataFormats/Candidate/interface/CompositeCandidate.h"
0017 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0018 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0019 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0020 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0021 #include <Math/VectorUtil.h>
0022
0023 class CAWZJetHelperUser {
0024 public:
0025 CAWZJetHelperUser(double massdropcut) : massdropcut_(massdropcut) {}
0026
0027 reco::CATopJetProperties operator()(reco::Jet const& ihardJet) const;
0028
0029 protected:
0030 double massdropcut_;
0031 };
0032
0033
0034
0035
0036
0037 class HLTCAWZTagFilter : public HLTFilter {
0038 public:
0039 explicit HLTCAWZTagFilter(const edm::ParameterSet&);
0040 ~HLTCAWZTagFilter() override;
0041 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0042 bool hltFilter(edm::Event&,
0043 const edm::EventSetup&,
0044 trigger::TriggerFilterObjectWithRefs& filterobject) const override;
0045
0046 private:
0047
0048
0049 edm::InputTag src_;
0050 edm::InputTag pfsrc_;
0051 const edm::EDGetTokenT<reco::BasicJetCollection> inputToken_;
0052 const edm::EDGetTokenT<reco::PFJetCollection> inputPFToken_;
0053 double minWMass_;
0054 double maxWMass_;
0055 double massdropcut_;
0056 };
0057
0058 inline reco::CATopJetProperties CAWZJetHelperUser::operator()(reco::Jet const& ihardJet) const {
0059 reco::CATopJetProperties properties;
0060
0061 reco::Jet::Constituents subjets = ihardJet.getJetConstituents();
0062 properties.nSubJets = subjets.size();
0063 properties.wMass = 999999.;
0064 properties.topMass = 999999.;
0065 properties.minMass = -1;
0066
0067 if (properties.nSubJets == 2) {
0068 sort(subjets.begin(), subjets.end(), [](auto const& t1, auto const& t2) { return t1->pt() > t2->pt(); });
0069
0070 reco::Jet::Constituent icandJet = subjets[0];
0071
0072 reco::Candidate::LorentzVector isubJet = icandJet->p4();
0073 double imass = isubJet.mass();
0074 double imw = ihardJet.mass();
0075
0076 if (imass / imw < massdropcut_) {
0077
0078 properties.wMass = imw;
0079 }
0080 }
0081
0082 return properties;
0083 }
0084
0085 #endif