rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2022_I2129461

tW single top-quark differential cross sections at 13 TeV
Experiment: CMS (LHC)
Inspire ID: 2129461
Status: VALIDATED
Authors:
  • Victor Rodriguez Bouza
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Measurement of inclusive and differential cross sections for single top quark production in association with a W boson in proton-proton collisions at $\sqrt{s}=$ 13 TeV

Measurements of the inclusive and normalised differential cross sections are presented for the production of single top quarks in association with a W boson in proton-proton collisions at a centre-of-mass energy of 13 TeV. The data used were recorded with the CMS detector at the LHC during 2016-2018, and correspond to an integrated luminosity of $138 fb^{-1}$. Events containing one electron and one muon in the final state are analysed. For the inclusive measurement, a multivariate discriminant, exploiting the kinematic properties of the events is used to separate the signal from the dominant $t\bar{t}$ background. A cross section of $79.2 \pm 0.9 (stat) ^{7.7}_{-8.0} (syst) \pm 1.2 (lumi) pb$ is obtained, consistent with the predictions of the standard model. For the differential measurements, a fiducial region is defined according to the detector acceptance, and the requirement of exactly one jet coming from the fragmentation of a bottom quark. The resulting distributions are unfolded to particle level and agree with the predictions at next-to-leading order in perturbative quantum chromodynamics.

Source code: CMS_2022_I2129461.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
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/ChargedLeptons.hh"
#include "Rivet/Projections/DressedLeptons.hh"
#include "Rivet/Projections/MissingMomentum.hh"
#include "Rivet/Projections/PromptFinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"

namespace Rivet {
  /// @brief Normalised differential cross sections for 13 TeV tW-channel single top-quark production
  class CMS_2022_I2129461 : public Analysis {
  public:

    /// Constructor
    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2022_I2129461);


    /// Book histograms and initialise projections before the run
    void init() {

      // final state of all stable particles
      Cut particle_cut = (Cuts::abseta < 5.0) and (Cuts::pT > 0.*MeV);
      FinalState fs(particle_cut);

      // select charged leptons
      ChargedLeptons charged_leptons(fs);

      // select final state photons for dressed lepton clustering
      IdentifiedFinalState photons(fs);
      photons.acceptIdPair(PID::PHOTON);

      // select final state prompt charged leptons
      PromptFinalState prompt_leptons(charged_leptons);
      prompt_leptons.acceptMuonDecays(true);
      prompt_leptons.acceptTauDecays(true);

      // select final state prompt photons
      PromptFinalState prompt_photons(photons);
      prompt_photons.acceptMuonDecays(true);
      prompt_photons.acceptTauDecays(true);

      // Dressed leptons from selected prompt charged leptons and photons
      Cut lepton_cut = ((Cuts::abseta < 2.4) and 
                        (Cuts::pT > 20.*GeV) and 
                        (((Cuts::abspid == PID::ELECTRON) and ((Cuts::abseta < 1.4442) or (Cuts::abseta > 1.566))) or (Cuts::abspid == PID::MUON)));

      DressedLeptons dressed_leptons(
        prompt_photons, prompt_leptons, 0.1,
        lepton_cut, true
      );
      declare(dressed_leptons, "DressedLeptons");

      // Jets
      VetoedFinalState fsForJets(fs);
      fsForJets.addVetoOnThisFinalState(dressed_leptons);
      declare(
        // excludes all neutrinos by default
        FastJets(fsForJets, FastJets::ANTIKT, 0.4),
        "Jets"
      );

      // pTmiss
      declare(MissingMomentum(fs), "MET");

      // book histograms
      book(_hist_norm_lep1_pt,            "d06-x01-y08");
      book(_hist_norm_lep1lep2jet1_pz,    "d07-x01-y08");
      book(_hist_norm_jet1_pt,            "d08-x01-y08");
      book(_hist_norm_lep1lep2jet1_m,     "d09-x01-y08");
      book(_hist_norm_lep1lep2_dphi,      "d10-x01-y08");
      book(_hist_norm_lep1lep2jet1met_mt, "d11-x01-y08");
    }


    /// @brief Perform the per-event analysis
    void analyze(const Event& event) {
      vector<DressedLepton> dressedLeptons = apply<DressedLeptons>(
        event,
        "DressedLeptons"
      ).dressedLeptons();

      // Require at least two dressed leptons
      if (dressedLeptons.size() < 2) vetoEvent;

      // Require that the leading ones are an electron and a muon of opposite charge
      if (abs(dressedLeptons.at(0).pid() + dressedLeptons.at(1).pid()) != 2) vetoEvent;

      // Require that the leading lepton has at least 25 GeV of pT
      if (dressedLeptons.at(0).pt() <= 25.*GeV) vetoEvent;

      // Require that all identified lepton pairs have at least 20 GeV of invariant mass
      for (int iL = 0; iL < int(dressedLeptons.size()); iL++) {
        for (int jL = iL + 1; jL < int(dressedLeptons.size()); jL++) {
          if ((dressedLeptons.at(iL).momentum() + dressedLeptons.at(jL).momentum()).mass() <= 20.*GeV) vetoEvent;
        }
      }

      // Jet object ID
      Cut jet_cut((Cuts::abseta < 2.4) and (Cuts::pT > 20.*GeV));
      vector<Jet> jets = apply<FastJets>(
        event,
        "Jets"
      ).jets(jet_cut);

      // ignore jets that overlap with dressed leptons within dR < 0.4 and collect them in loose jets and jets categories
      Jets cleanedLooseJets;
      Jets cleanedJets;
      idiscardIfAnyDeltaRLess(jets, dressedLeptons, 0.4);
      for (const Jet& jet: jets) {
        if (jet.pt() > 30.*GeV) cleanedJets.push_back(jet);
        else                    cleanedLooseJets.push_back(jet);
      }

      if (cleanedJets.size() != 1)         vetoEvent; // select events with exactly one jet...
      if (not cleanedJets.at(0).bTagged()) vetoEvent; // ...that must be b-tagged...
      if (cleanedLooseJets.size() > 0)     vetoEvent; // ...and no loose jets

      // fill the histograms
      _hist_norm_lep1_pt->fill(min(max(dressedLeptons.at(0).pt(), 26.*GeV), 149.*GeV) / GeV);
      FourMomentum llj = dressedLeptons.at(0).momentum() + dressedLeptons.at(1).momentum() + cleanedJets.at(0).momentum();
      _hist_norm_lep1lep2jet1_pz->fill(min(max(abs((llj).pz()), 1.*GeV), 449.*GeV) / GeV);
      _hist_norm_jet1_pt->fill(min(max(cleanedJets.at(0).pt(), 31.*GeV), 149.*GeV) / GeV);
      _hist_norm_lep1lep2jet1_m->fill(min(max((llj).mass(), 51.*GeV), 399.*GeV) / GeV);
      _hist_norm_lep1lep2_dphi->fill( abs(deltaPhi(dressedLeptons.at(0).phi(), dressedLeptons.at(1).phi())) / PI );

      Vector3 tmpPtmiss3  = apply<MissingMomentum>(event, "MET").vectorPt();
      tmpPtmiss3.setZ(0.);
      FourMomentum ptmiss = FourMomentum(tmpPtmiss3.mod(), -tmpPtmiss3.x(), -tmpPtmiss3.y(), 0.);
      FourMomentum tmpM   = llj + ptmiss;
      _hist_norm_lep1lep2jet1met_mt->fill(min(max(sqrt(tmpM.E2() - tmpM.pz2()), 101.*GeV), 499.*GeV) / GeV);
    }

    /// @brief Normalise histograms after the run
    void finalize() {
      normalize(_hist_norm_lep1_pt);
      normalize(_hist_norm_lep1lep2jet1_pz);
      normalize(_hist_norm_jet1_pt);
      normalize(_hist_norm_lep1lep2jet1_m);
      normalize(_hist_norm_lep1lep2_dphi);
      normalize(_hist_norm_lep1lep2jet1met_mt);
    }

  private:
    // Declaration of histograms
    Histo1DPtr _hist_norm_lep1_pt;
    Histo1DPtr _hist_norm_lep1lep2jet1_pz;
    Histo1DPtr _hist_norm_jet1_pt;
    Histo1DPtr _hist_norm_lep1lep2jet1_m;
    Histo1DPtr _hist_norm_lep1lep2_dphi;
    Histo1DPtr _hist_norm_lep1lep2jet1met_mt;
  };

  RIVET_DECLARE_PLUGIN(CMS_2022_I2129461);
}