Rivet Analyses Reference

CMS_2016_I1454211

Boosted $t\bar{t}$ in $pp$ collisions at $\sqrt{s} = 8 \text{TeV}$
Experiment: CMS (LHC)
Inspire ID: 1454211
Status: VALIDATED
Authors:

Salvatore Rappoccio
Maral Alyari
Julia Thom
Louise Skinnari
Susan Dittmer
Matthew Bellis

References:

CMS-PAS-TOP-14-012 (PAS)
CERN-EP-2016-078 (paper)

Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:

$pp$ QCD interactions at $\sqrt{s} = 8 \text{TeV} with $t\bar{t}$ process. Data collected by CMS during the year 2012. Boosted topology restricts phase space, requiring high statistics in a single run or custom merging of the YODA files.

The cross section for pair production of top quarks with high transverse momenta ($p_\mathrm{T} > 400$ \text{GeV}) is measured in 19.7 fb$^{-1}$ of $\mathrm{pp}$ collisions, collected with the CMS detector at $\sqrt{s} = 8 \text{TeV}$. The measurement is performed for lepton+jets events, where one top quark decays according to $t \rightarrow Wb \rightarrow \ell \nu b$, with $\ell$ denoting an electron or muon, and the second top quark decays to an hadronic final state and is reconstructed as a single large-radius jet and identified as a top quark candidate using jet substructure techniques. Integrated cross sections, as well as differential cross sections as a function of the top quark $p_\mathrm{T}$ and rapidity, are measured both at particle level within a fiducial region resembling the detector-level selections and at parton level. RIVET: This analysis is to be run on $t\bar{t}$ Monte Carlo. It utilizes the PartonicTops projection, which assumes top quarks in the event record. The analysis has been validated with Powheg+Pythia6. The parton-level phase space is defined by requiring two PartonicTops. Exactly one PartonicTop must decay directly to a muon or electron (no intermediate tau), and exactly one PartonicTop decays hadronically. For $t\bar{t}$ Monte Carlo, this is equivalent to requiring the event to be semileptonic at parton level. The parton-level top quark is defined as the hadronically decaying top. The parton-level top quark is required to have $p_\mathrm{T} > 400 \text{GeV}$. The particle-level phase space is defined using the following object definitions: - Lepton: A dressed electron or muon, meaning the lepton has been clustered with all photons within a cone of $R=0.1$. The DressedLepton projection is used to construct the dressed lepton. The lepton is required to have $p_\mathrm{T} > 45$ GeV and $|\eta| < 2.1$. - B Jet Candidate: Gen AK5 jets are formed by clustering the final state particles in the event using the anti-$k_{T}$ algorithm with distance parameter $R=0.5$. Neutrinos are excluded from the clustering, as are any particles included in the dressed lepton. The gen AK5 jet is required to have $p_{T} > 30$ GeV and $|\eta| < 2.4$. Gen AK5 jets in the same hemisphere as the lepton ($\Delta R(e/\mu, \mathrm{jet}) < \pi/2$) are defined as $b$-jet candidates. - Top Jet Candidate: Gen CA8 jets are formed by clustering the final state particles in the event using the Cambridge-Aachen algorithm with distance parameter $R=0.8$. Neutrinos are excluded from the clustering, as are any particles included in the dressed lepton. The gen CA8 jet is required to have $p_{T} > 30$ GeV and $|\eta| < 2.4$. Gen CA8 jets which have $p_{T} > 400$ GeV, 140 GeV $<$ mass $<$ 250 GeV, and are in the opposite hemisphere from the lepton ($\Delta{\rm R(e/\mu, jet)} > \pi/2$) are defined as top jet candidates. The particle-level phase space is defined by requiring $\geq 1$ b jet candidate, $\geq 1$ top jet candidate, and exactly one lepton. This is in addition to the parton-level semileptonic requirement. The highest-$p_\mathrm{T}$ top jet candidate is defined as the particle-level $t$ jet.

Source code: CMS_2016_I1454211.cc

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/PartonicTops.hh"
#include "Rivet/Projections/DressedLeptons.hh"
#include "Rivet/Projections/IdentifiedFinalState.hh"
#include "Rivet/Projections/PromptFinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/InvMassFinalState.hh"
#include "Rivet/Projections/MissingMomentum.hh"

namespace Rivet {


  /// Boosted ttbar in pp collisions at sqrtS = 8 TeV
  /// @todo Use persistent weight counters
  class CMS_2016_I1454211 : public Analysis {
  public:
    
    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1454211);
    
    
    // Set up projections and book histograms
    void init() {
      
      // Get options particle-level only.
      _mode = 0;
      if ( getOption("TMODE") == "PARTICLE" ) _mode = 0;
      if ( getOption("TMODE") == "BOTH" ) _mode = 1;
      
      // Complete final state
      FinalState fs;
      
      // Partonic tops
      // Need these for flavour determination, even if only plotting particle-level
      declare(PartonicTops(PartonicTops::DecayMode::ELECTRON, false), "ElectronPartonTops");
      declare(PartonicTops(PartonicTops::DecayMode::MUON, false),     "MuonPartonTops");
      declare(PartonicTops(PartonicTops::DecayMode::HADRONIC),        "HadronicPartonTops");
      
      // Projection for electrons and muons
      IdentifiedFinalState photons(fs, PID::PHOTON);
      
      const Cut leptonCuts = Cuts::pt > 45*GeV && Cuts::abseta < 2.1;
      
      IdentifiedFinalState el_id(fs, {{PID::ELECTRON, -PID::ELECTRON}});
      PromptFinalState electrons(el_id);
      DressedLeptons dressed_electrons(photons, electrons, 0.1, leptonCuts);
      declare(dressed_electrons, "DressedElectrons");
      
      IdentifiedFinalState mu_id(fs, {{PID::MUON, -PID::MUON}});
      PromptFinalState muons(mu_id);
      DressedLeptons dressed_muons(photons, muons, 0.1, leptonCuts);
      declare(dressed_muons, "DressedMuons");

      // Projection for jets
      VetoedFinalState fs_jets(fs);
      fs_jets.addVetoOnThisFinalState(dressed_muons);
      fs_jets.addVetoOnThisFinalState(dressed_electrons);
      fs_jets.vetoNeutrinos();
      declare(FastJets(fs_jets, FastJets::ANTIKT, 0.5), "ak5jets");
      declare(FastJets(fs_jets, FastJets::CAM, 0.8), "ca8jets");

      if (_mode == 1) {
        book(_hEl_topPt_parton          , "d01-x01-y01"); // dsigma/dpt(top quark), el ch
        book(_hEl_topY_parton           , "d03-x01-y01"); // dsigma/dy(top quark), el ch
        book(_hMu_topPt_parton          , "d05-x01-y01"); // dsigma/dpt(top quark), mu ch
        book(_hMu_topY_parton           , "d07-x01-y01"); // dsigma/dy(top quark), mu ch
        book(_hComb_topPt_parton        , "d09-x01-y01"); // dsigma/dpt(top quark), comb ch
        book(_hComb_topY_parton         , "d11-x01-y01"); // dsigma/dy(top quark), comb ch

        book(_hEl_topPt_parton_norm     , "d13-x01-y01"); // 1/sigma dsigma/dpt(top quark), el ch
        book(_hEl_topY_parton_norm      , "d15-x01-y01"); // 1/sigma dsigma/dy(top quark), el ch
        book(_hMu_topPt_parton_norm     , "d17-x01-y01"); // 1/sigma dsigma/dpt(top quark), mu ch
        book(_hMu_topY_parton_norm      , "d19-x01-y01"); // 1/sigma dsigma/dy(top quark), mu ch
        book(_hComb_topPt_parton_norm   , "d21-x01-y01"); // 1/sigma dsigma/dpt(top quark), comb ch
        book(_hComb_topY_parton_norm    , "d23-x01-y01"); // 1/sigma dsigma/dy(top quark), comb ch
      }

      book(_hEl_topPt_particle        , "d02-x01-y01"); // dsigma/dpt(top jet), el ch
      book(_hEl_topY_particle         , "d04-x01-y01"); // dsigma/dy(top jet), el ch
      book(_hMu_topPt_particle        , "d06-x01-y01"); // dsigma/dpt(top jet), mu ch
      book(_hMu_topY_particle         , "d08-x01-y01"); // dsigma/dy(top jet), mu ch
      book(_hComb_topPt_particle      , "d10-x01-y01"); // dsigma/dpt(top jet), comb ch
      book(_hComb_topY_particle       , "d12-x01-y01"); // dsigma/dy(top jet), comb ch
      book(_hEl_topY_particle_norm    , "d16-x01-y01"); // 1/sigma dsigma/dy(top jet), el ch
      book(_hEl_topPt_particle_norm   , "d14-x01-y01"); // 1/sigma dsigma/dpt(top jet), el ch
      book(_hComb_topY_particle_norm  , "d24-x01-y01"); // 1/sigma dsigma/dy(top jet), comb ch
      book(_hMu_topPt_particle_norm   , "d18-x01-y01"); // 1/sigma dsigma/dpt(top jet), mu ch
      book(_hMu_topY_particle_norm    , "d20-x01-y01"); // 1/sigma dsigma/dy(top jet), mu ch
      book(_hComb_topPt_particle_norm , "d22-x01-y01"); // 1/sigma dsigma/dpt(top jet), comb ch

      book(_hMu_cutflow , "mu_cutflow", 7, -0.5, 6.5);
      book(_hEl_cutflow , "el_cutflow", 7, -0.5, 6.5);
    }


    // per event analysis
    void analyze(const Event& event) {

       // Total-events cutflow entries
      _hMu_cutflow->fill(0.);
      _hEl_cutflow->fill(0.);

      // Do parton-level selection and channel determination
      // Note that channel determination relies on partonic info, even for the particle-level tops
      int partonCh = 0; //0 non-semi-lep, 1 muon, 2 electron
      const Particles muonpartontops = apply<ParticleFinder>(event, "MuonPartonTops").particlesByPt();
      const Particles electronpartontops = apply<ParticleFinder>(event, "ElectronPartonTops").particlesByPt();
      if (electronpartontops.size() == 0 && muonpartontops.size() == 1) partonCh = 1;
      else if (electronpartontops.size() == 1 && muonpartontops.size() == 0) partonCh = 2;
      else vetoEvent;
      const Particles hadronicpartontops = apply<ParticleFinder>(event, "HadronicPartonTops").particlesByPt();
      if (hadronicpartontops.size() != 1) vetoEvent;

      if (partonCh == 1) _hMu_cutflow->fill(1.); // muon at parton level
      if (partonCh == 2) _hEl_cutflow->fill(1.); // electron at parton level

      // Get hadronic parton-level top
      const FourMomentum& partonTopP4 = hadronicpartontops.front();

      // Do particle-level selection and channel determination
      const DressedLeptons& dressed_electrons = apply<DressedLeptons>(event, "DressedElectrons");
      const DressedLeptons& dressed_muons = apply<DressedLeptons>(event, "DressedMuons");

      bool passParticleLep = false, passParticleTop = false;
      FourMomentum lepton, particleTopP4;

      if (partonCh == 1 && dressed_muons.dressedLeptons().size() == 1 && dressed_electrons.dressedLeptons().size() == 0) {
        passParticleLep = true;
        _hMu_cutflow->fill(3.); //muon at particle level
        lepton = dressed_muons.dressedLeptons()[0].momentum();
      }
      if (partonCh == 2 && dressed_muons.dressedLeptons().size() == 0 && dressed_electrons.dressedLeptons().size() == 1) {
        passParticleLep = true;
        _hEl_cutflow->fill(3.); //electron at particle level
        lepton = dressed_electrons.dressedLeptons()[0].momentum();
      }

      if (passParticleLep) {

        // Jet cuts
        Cut jetCuts = Cuts::pt > 30*GeV && Cuts::abseta < 2.4;
        Jets genBjets, genTjets;
        int nGenBjets = 0, nGenTjets = 0;
        
        const FastJets& AK5jets = apply<FastJets>(event, "ak5jets");
        for (const Jet& jet : AK5jets.jetsByPt(jetCuts)) {
          if (deltaR(jet, lepton) > M_PI / 2.0) continue;
          if (deltaR(jet, lepton) < 0.1) continue;
          genBjets.push_back(jet);
          nGenBjets += 1;
        }
        
        const FastJets& CA8jets = apply<FastJets>(event, "ca8jets");
        for (const Jet& jet : CA8jets.jetsByPt(jetCuts)) {
          if (deltaR(jet, lepton) < M_PI / 2.0) continue;
          if (jet.mass() < 140*GeV) continue;
          if (jet.mass() > 250*GeV) continue;
          genTjets.push_back(jet);
          nGenTjets += 1;
        }
        
        if (nGenBjets >=1) {
          if (_mode == 1) {
            if (partonCh == 1) _hMu_cutflow->fill(4.); // muon at parton level
            if (partonCh == 2) _hEl_cutflow->fill(4.); // electron at parton level
          }
          if (nGenTjets >= 1) {
            passParticleTop = true;
            if (_mode == 1) {
              if (partonCh == 1) _hMu_cutflow->fill(5.); // muon at parton level
              if (partonCh == 2) _hEl_cutflow->fill(5.); // electron at parton level
            }
            particleTopP4 = genTjets[0];
          }
        }
      }
      
      const double weight = 1.0;        
      if (partonCh == 1) {
        _nMu += weight;

        if (_mode == 1) {
          // protect against unphysical partons
          if (partonTopP4.E() < 0) {
            MSG_WARNING("Top parton with negative energy! Vetoing event. Try turning off partonic tops?");
            vetoEvent;
          }
          
          _hMu_topPt_parton->fill(partonTopP4.pT()/GeV, weight);
          _hMu_topPt_parton_norm->fill(partonTopP4.pT()/GeV, weight);
          _hComb_topPt_parton->fill(partonTopP4.pT()/GeV, weight);
          _hComb_topPt_parton_norm->fill(partonTopP4.pT()/GeV, weight);
          
          if (partonTopP4.pT() >= 400*GeV) {
            _nPassParton_mu += weight;
            _hMu_cutflow->fill(2.);
            _hMu_topY_parton->fill(partonTopP4.rapidity(), weight);
            _hMu_topY_parton_norm->fill(partonTopP4.rapidity(), weight);
            _hComb_topY_parton->fill(partonTopP4.rapidity(), weight);
            _hComb_topY_parton_norm->fill(partonTopP4.rapidity(), weight);
          }
        }
         
        if (passParticleTop) {
          _hMu_topPt_particle->fill(particleTopP4.pT()/GeV, weight);
          _hMu_topPt_particle_norm->fill(particleTopP4.pT()/GeV, weight);
          _hComb_topPt_particle->fill(particleTopP4.pT()/GeV, weight);
          _hComb_topPt_particle_norm->fill(particleTopP4.pT()/GeV, weight);
          
          if (particleTopP4.pT() >= 400*GeV) {
            _nPassParticle_mu += weight;
            _hMu_cutflow->fill(6.);
            _hMu_topY_particle->fill(particleTopP4.rapidity(), weight);
            _hMu_topY_particle_norm->fill(particleTopP4.rapidity(), weight);
            _hComb_topY_particle->fill(particleTopP4.rapidity(), weight);
            _hComb_topY_particle_norm->fill(particleTopP4.rapidity(), weight);
          }
        }
      }

      if (partonCh == 2){
        _nEl += weight;
        if (_mode == 1) {
          _hEl_topPt_parton->fill(partonTopP4.pT()/GeV, weight);
          _hEl_topPt_parton_norm->fill(partonTopP4.pT()/GeV, weight);
          _hComb_topPt_parton->fill(partonTopP4.pT()/GeV, weight);
          _hComb_topPt_parton_norm->fill(partonTopP4.pT()/GeV, weight);
          
          if (partonTopP4.pT() >= 400*GeV) {
            _nPassParton_el += weight;
            _hEl_cutflow->fill(2.);
            _hEl_topY_parton->fill(partonTopP4.rapidity(), weight);
            _hEl_topY_parton_norm->fill(partonTopP4.rapidity(), weight);
            _hComb_topY_parton->fill(partonTopP4.rapidity(), weight);
            _hComb_topY_parton_norm->fill(partonTopP4.rapidity(), weight);
          }
        }
      
        if (passParticleTop) {
          _hEl_topPt_particle->fill(particleTopP4.pT()/GeV, weight);
          _hEl_topPt_particle_norm->fill(particleTopP4.pT()/GeV, weight);
          _hComb_topPt_particle->fill(particleTopP4.pT()/GeV, weight);
          _hComb_topPt_particle_norm->fill(particleTopP4.pT()/GeV, weight);
          
          if (particleTopP4.pT() >= 400*GeV) {
            _nPassParticle_el += weight;
            _hEl_cutflow->fill(6.);
            _hEl_topY_particle->fill(particleTopP4.rapidity(), weight);
            _hEl_topY_particle_norm->fill(particleTopP4.rapidity(), weight);
            _hComb_topY_particle->fill(particleTopP4.rapidity(), weight);
            _hComb_topY_particle_norm->fill(particleTopP4.rapidity(), weight);
          }
        }
      }
    }
    
    void finalize() {

      normalize(_hMu_topPt_particle_norm); normalize(_hMu_topY_particle_norm); normalize(_hEl_topPt_particle_norm);
      normalize(_hEl_topY_particle_norm); normalize(_hComb_topPt_particle_norm); normalize(_hComb_topY_particle_norm, 1.0, false);
      
      const double sf = crossSection() / femtobarn / sumOfWeights();
      scale(_hMu_topPt_particle, sf);
      scale(_hEl_topPt_particle, sf);
      scale(_hMu_topY_particle, sf);
      scale(_hEl_topY_particle, sf);
      scale(_hComb_topPt_particle, sf);
      scale(_hComb_topY_particle, sf);
      
      if (_mode == 1) {
        normalize(_hMu_topPt_parton_norm); normalize(_hMu_topY_parton_norm); normalize(_hEl_topPt_parton_norm);
        normalize(_hEl_topY_parton_norm); normalize(_hComb_topPt_parton_norm); normalize(_hComb_topY_parton_norm, 1.0, false);
        scale(_hMu_topPt_parton, sf);
        scale(_hEl_topPt_parton, sf);
        scale(_hMu_topY_parton, sf);
        scale(_hEl_topY_parton, sf);
        scale(_hComb_topPt_parton, sf);
        scale(_hComb_topY_parton, sf);
      }
    }

  protected:

    size_t _mode;

  private:

    Histo1DPtr _hMu_topPt_parton, _hMu_topY_parton, _hEl_topPt_parton, _hEl_topY_parton, _hComb_topPt_parton, _hComb_topY_parton;
    Histo1DPtr _hMu_topPt_particle, _hMu_topY_particle, _hEl_topPt_particle, _hEl_topY_particle, _hComb_topPt_particle, _hComb_topY_particle;
    Histo1DPtr _hMu_topPt_parton_norm, _hMu_topY_parton_norm, _hEl_topPt_parton_norm, _hEl_topY_parton_norm, _hComb_topPt_parton_norm, _hComb_topY_parton_norm;
    Histo1DPtr _hMu_topPt_particle_norm, _hMu_topY_particle_norm, _hEl_topPt_particle_norm, _hEl_topY_particle_norm, _hComb_topPt_particle_norm, _hComb_topY_particle_norm;
    Histo1DPtr _hMu_cutflow, _hEl_cutflow;

    double _nMu = 0., _nEl = 0.;
    double _nPassParton_mu = 0.,_nPassParton_el = 0.;
    double _nPassParticle_mu = 0., _nPassParticle_el = 0.;

  };


  // The hook for the plugin system
  RIVET_DECLARE_PLUGIN(CMS_2016_I1454211);

}

Edit this page on GitLab