Rivet/Tools/RivetHTT.hh

0001 // -*- C++ -*-
0002 #ifndef RIVET_RIVETHTT_HH
0003 #define RIVET_RIVETHTT_HH
0004
0005 #include "Rivet/Jet.hh"
0006 #include "HEPTopTagger/HEPTopTagger.hh"
0007
0008 namespace Rivet {
0009
0010
0011   /// HTT operating mode
0012   enum class HTTMode {
0013     EARLY_MASSRATIO_SORT_MASS,     // applies 2D mass plane requirements then select the candidate which minimizes |m_cand-mt|
0014     LATE_MASSRATIO_SORT_MASS,      // selects the candidate which minimizes |m_cand-mt|
0015     EARLY_MASSRATIO_SORT_MODDJADE, // applies the 2D mass plane requirements then select the candidate with highest jade distance
0016     LATE_MASSRATIO_SORT_MODDJADE,  // selects the candidate with highest modified jade distance
0017     TWO_STEP_FILTER                // only analyzes the candidate built with the highest pT(t) after unclustering
0018   };
0019
0020
0021   /// Wrapper class for configuring use of HEPTopTagger
0022   class HTT {
0023   public:
0024
0025     struct InputParameters {
0026
0027       /// HTT execution mode
0028       Mode mode = HTTMode::EARLY_MASSRATIO_SORT_MASS;
0029
0030       /// @name Optimal-R parameters
0031       /// @{
0032       /// Use the optimal-R finder (default), or use fixed R
0033       bool do_optimalR = true;
0034       double optimalR_min = 0.5; // min jet size
0035       double optimalR_step = 0.1; // step size
0036       double optimalR_threshold = 0.2; // step size
0037       /// @}
0038
0039       /// @name Mass-drop declustering
0040       /// @{
0041       double mass_drop = 0.8;
0042       double max_subjet_mass = 30*GeV;
0043       /// @}
0044
0045       /// @name Filtering
0046       /// @{
0047       unsigned int filt_N = 5; // set_nfilt
0048       double filtering_R = 0.3; // max subjet distance for filtering
0049       double filtering_minpt = 0.; // min subjet pt for filtering
0050       /// @}
0051
0052       /// Jet algorithm used for filtering
0053       JetAlg filtering_algorithm = JetAlg::CA;
0054
0055       /// Reclustering jet algorithm
0056       JetAlg reclustering_algorithm = JetAlg::CA;
0057
0058       /// @name Top-mass ranges
0059       /// @{
0060       /// @todo Take from a central set of constants
0061       double top_mass = 172.3*GeV;
0062       /// @todo Take from a central set of constants
0063       double W_mass = 80.4*GeV;
0064       double Mtop_min = 150*GeV;
0065       double Mtop_max = 200*GeV; //set_top_range(min,max)
0066       /// @}
0067
0068       /// @name Top-mass ratio range
0069       /// @{
0070       double fw = 0.15;
0071       double mass_ratio_range_min = (1.-fw)*W_mass/top_mass;
0072       double mass_ratio_range_max = (1.+fw)*W_mass/top_mass;
0073       /// @}
0074
0075       /// @name Mass-ratio cuts
0076       /// @{
0077       double m23cut = 0.35;
0078       double m13cutmin = 0.2;
0079       double m13cutmax = 1.3;
0080       /// @}
0081
0082       /// @name Pruning cuts
0083       /// @{
0084       double prune_zcut = 0.1;
0085       double prune_rcut = 0.5;
0086       /// @}
0087
0088     };
0089
0090
0091
0092     /// Constructor without arguments
0093     HTT() {}
0094
0095     /// Constructor with arguments
0096     HTT(HTT::InputParameters& params) {
0097       setParams(params);
0098     }
0099
0100     // /// Destructor
0101     // ~HTT() {}
0102
0103     /// Set the tagging parameters
0104     void setParams(HTT::InputParameters& params);
0105
0106     /// Run the top tagger on a given jet
0107     void calc(Jet& jet);
0108
0109     /// Top jet
0110     const Jet topJet() const;
0111     /// The bottom jet inside the top
0112     const Jet bJet() const;
0113     /// The W jet inside the top
0114     const Jet wJet() const;
0115     /// Leading subjet from W
0116     const Jet w1Jet() const;
0117     /// Second leading subjet from W
0118     const Jet w2Jet() const;
0119
0120     /// Top jet, as a pseudojet
0121     const PseudoJet& topJet() const;
0122     /// The bottom jet inside the top, as a pseudojet
0123     const PseudoJet& bJet() const;
0124     /// The W jet inside the top, as a pseudojet
0125     const PseudoJet& wJet() const;
0126     /// Leading subjet from W, as a pseudojet
0127     const PseudoJet& w1Jet() const;
0128     /// Second leading subjet from W, as a pseudojet
0129     const PseudoJet& w2Jet() const;
0130
0131     /// pT-ordered subjets
0132     const Jets& subjets() const;
0133
0134     // /// Print tagger information
0135     // void info() const;
0136     // /// Print tagger settings
0137     // void settings() const;
0138
0139     /// The pruned mass
0140     double prunedMass() const;
0141
0142     // The unfiltered mass
0143     double unfilteredMass() const;
0144
0145     /// Difference between the reco top mass and the true top mass
0146     double deltaTopMass() const;
0147
0148     /// Is the jet tagged?
0149     bool isTopTagged() const;
0150
0151     /// Was the top-mass window requirement passed?
0152     bool passedMassCutTop() const;
0153
0154     /// 2D mass plane requirements passed?
0155     bool passedMassCut2D() const;
0156
0157
0158   private:
0159
0160     fastjet::HEPTopTagger::HEPTopTagger _tagger;
0161
0162   };
0163
0164
0165   /// Below can be moved to HTT.cc at some point
0166
0167   void HTT::setParams(HTT::InputParameters& params) {
0168     _tagger = fastjet::HEPTopTagger::HEPTopTagger();
0169
0170     // Optimal R
0171     _tagger.do_optimalR(params.do_optimalR);
0172     _tagger.set_optimalR_min(params.optimalR_min);
0173     _tagger.set_optimalR_step(params.optimalR_step);
0174     _tagger.set_optimalR_threshold(params.optimalR_threshold);
0175
0176     // Candidate selection
0177     fastjet::HEPTopTagger::Mode mode;
0178     if (params.mode == HTTMode::EARLY_MASSRATIO_SORT_MASS) {
0179       mode = fastjet::HEPTopTagger::EARLY_MASSRATIO_SORT_MASS;
0180     } else if (params.mode == HTTMode::LATE_MASSRATIO_SORT_MASS) {
0181       mode = fastjet::HEPTopTagger::LATE_MASSRATIO_SORT_MASS;
0182     } else if (params.mode == HTTMode::EARLY_MASSRATIO_SORT_MODDJADE) {
0183       mode = fastjet::HEPTopTagger::EARLY_MASSRATIO_SORT_MODDJADE;
0184     } else if (params.mode == HTTMode::LATE_MASSRATIO_SORT_MODDJADE) {
0185       mode = fastjet::HEPTopTagger::LATE_MASSRATIO_SORT_MODDJADE;
0186     } else {
0187       mode = fastjet::HEPTopTagger::TWO_STEP_FILTER;
0188     }
0189     _tagger.set_mode(mode);
0190     _tagger.set_mt(params.top_mass);
0191     _tagger.set_mw(params.W_mass);
0192     _tagger.set_top_mass_range(params.Mtop_min, params.Mtop_max);
0193     _tagger.set_fw(params.fw);
0194     _tagger.set_mass_ratio_range(params.mass_ratio_range_min, params.mass_ratio_range_max);
0195     _tagger.set_mass_ratio_cut(params.m23cut, params.m13cutmin, params.m13cutmax);
0196
0197     // Filtering
0198     _tagger.set_filtering_n(params.filt_N);
0199     _tagger.set_filtering_R(params.filtering_R);
0200     _tagger.set_filtering_minpt_subjet(params.filtering_minpt);
0201
0202     fastjet::JetAlgorithm algo = fastjet::antikt_algorithm;
0203     if (params.filtering_algorithm == Algo::CA) algo = fastjet::cambridge_algorithm;
0204     else if (params.filtering_algorithm == Algo::KT) algo = fastjet::kt_algorithm;
0205     _tagger.set_filtering_jetalgorithm(algo);
0206
0207     // Reclustering
0208     algo = fastjet::antikt_algorithm;
0209     if (params.reclustering_algorithm == Algo::CA) algo = fastjet::cambridge_algorithm;
0210     else if (params.reclustering_algorithm == Algo::KT) algo = fastjet::kt_algorithm;
0211     _tagger.set_reclustering_jetalgorithm(algo);
0212
0213     // Mass-drop
0214     _tagger.set_mass_drop_threshold(params.mass_drop);
0215     _tagger.set_mass_drop_threshold(params.mass_drop);
0216
0217     // Pruning
0218     _tagger.set_pruning_rcut_factor(params.prune_rcut);
0219     _tagger.set_pruning_zcut(params.prune_zcut);
0220   }
0221
0222
0223   void HTT::calc(Jet& jet) {
0224     _tagger.run(jet.pseudojet());
0225   }
0226
0227
0228   const Jet HTT::topJet() const { return Jet(topPjet()); }
0229   const Jet HTT::bJet() const { return Jet(bPjet()); }
0230   const Jet HTT::wJet() const { return Jet(wPjet()); }
0231   const Jet HTT::w1Jet() const { return Jet(w1Pjet()); }
0232   const Jet HTT::w2Jet() const { return Jet(w2Pjet()); }
0233
0234
0235   const PseudoJet& HTT::topPjet() const { return _tagger.t(); }
0236   const PseudoJet& HTT::bPjet() const { return _tagger.b(); }
0237   const PseudoJet& HTT::wPjet() const { return _tagger.W(); }
0238   const PseudoJet& HTT::w1Pjet() const { return _tagger.W1(); }
0239   const PseudoJet& HTT::w2Pjet() const { return _tagger.W2(); }
0240
0241
0242   Jets HTT::subjets() const {
0243     Jets rtn;
0244     rtn.reserve(3);
0245     rtn.emplace_back(_tagger.j1());
0246     rtn.emplace_back(_tagger.j2());
0247     rtn.emplace_back(_tagger.j3());
0248     return rtn;
0249   }
0250
0251
0252   // void HTT::info() const { _tagger.get_info(); }
0253   // void HTT::settings() const { _tagger.get_setting(); }
0254
0255   double HTT::prunedMass() const { return _tagger.pruned_mass(); }
0256
0257   double HTT::unfilteredMass() const { return _tagger.unfiltered_mass(); }
0258
0259   double HTT::deltaTopMass() const { return _tagger.delta_top(); }
0260
0261   bool HTT::isTopTagged() const { return _tagger.is_tagged(); }
0262
0263   bool HTT::passedMassCutTop() const { return _tagger.is_maybe_top(); }
0264
0265   bool HTT::passedMassCut2D() const { return _tagger.is_masscut_passed(); }
0266
0267   /// Above can be moved to RivetHTT.cc at some point
0268
0269
0270 }
0271
0272 #endif