rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2017_I1634835

Measurement of associated Z + charm production in proton-proton collisions at sqrts = 8 TeV
Experiment: CMS (LHC)
Inspire ID: 1634835
Status: VALIDATED
Authors:
  • Hannes Jung
  • Juan Pablo Fernandez
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
    No run details listed

A study of the associated production of a Z boson and a charm quark jet (Z+c) and a comparison to production with a b quark jet (Z+b), in pp collisions at a centre-of-mass energy of 8 TeV are presented. The analysis uses a data sample cor- responding to an integrated luminosity of 19.7 fb$^{-1}$, collected with the CMS detector at the CERN LHC. The Z boson candidates are identified through their decays into pairs of electrons or muons. Jets originating from heavy flavour quarks are iden- tified using semileptonic decays of c or b flavoured hadrons and hadronic decays of charm hadrons. The measurements are performed in the kinematic region with two leptons with $p_T^l > 20$ GeV, $|\eta^l| < 2.1$, $71 < m_{ll} < 111$ GeV, and heavy flavour jets with $p_T^{jet} > 25$ GeV and $|\eta^{jet}| <2.5$. The Z + c production cross section and the cross section ratio are also measured as a function of the transverse momentum of the Z boson and of the heavy flavour jet.

Source code: CMS_2017_I1634835.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
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/ZFinder.hh"
#include "Rivet/Tools/Cuts.hh"

namespace Rivet {

/// @brief Z+charm at 8 TeV
class CMS_2017_I1634835 : public Analysis {
   public:
    /// Constructor
    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1634835);

    /// @name Analysis methods
    ///@{

    /// Book histograms and initialise projections before the run
    void init() {
        // Initialise and register projections
        FinalState fs;

        ZFinder zeeFinder(fs, Cuts::abseta < 2.1 && Cuts::pT > 20 * GeV, PID::ELECTRON, 71.0 * GeV, 111.0 * GeV, 0.1);
        declare(zeeFinder, "ZeeFinder");

        ZFinder zmumuFinder(fs, Cuts::abseta < 2.1 && Cuts::pT > 20 * GeV, PID::MUON, 71.0 * GeV, 111.0 * GeV, 0.1);
        declare(zmumuFinder, "ZmumuFinder");

        VetoedFinalState jetConstits(fs);
        jetConstits.addVetoOnThisFinalState(zeeFinder);
        jetConstits.addVetoOnThisFinalState(zmumuFinder);

        FastJets akt04Jets(jetConstits, FastJets::ANTIKT, 0.4);
        declare(akt04Jets, "AntiKt04Jets");

        book(_h_Z_pt_cjet, 4, 1, 1);
        book(_h_pt_cjet, 5, 1, 1);
        book(_h_Z_pt_bjet, "TMP/Z_pt_bjet", refData(4, 1, 2));
        book(_h_pt_bjet, "TMP/pt_bjet", refData(5, 1, 2));
        // book ratio histos

        book(_h_R_Z_pt_cjet, 4, 1, 2);
        book(_h_R_jet_pt_cjet, 5, 1, 2);

        book(counter_ee, "TMP/counter_ee");
        book(counter_mm, "TMP/counter_mm");
    }

    /// Perform the per-event analysis
    void analyze(const Event& event) {
        const ZFinder& zeeFS = apply<ZFinder>(event, "ZeeFinder");
        const ZFinder& zmumuFS = apply<ZFinder>(event, "ZmumuFinder");

        const Particles& zees = zeeFS.bosons();
        const Particles& zmumus = zmumuFS.bosons();

        // We did not find exactly one Z. No good.
        if (zees.size() + zmumus.size() != 1) {
            MSG_DEBUG("Did not find exactly one good Z candidate");
            vetoEvent;
        }

        Particles leptons;
        Particle zcand;
        if (zees.size() == 1) {
            leptons = zeeFS.constituentLeptons();
            zcand = zees[0];
            counter_ee->fill();
        }
        if (zmumus.size() == 1) {
            leptons = zmumuFS.constituentLeptons();
            zcand = zmumus[0];
            counter_mm->fill();
        }

        // Cluster jets
        // NB. Veto has already been applied on leptons and photons used for dressing
        const FastJets& fj = apply<FastJets>(event, "AntiKt04Jets");
        Jets goodjets = fj.jetsByPt(Cuts::abseta < 2.5 && Cuts::pT > 25*GeV);
        idiscardIfAnyDeltaRLess(goodjets, leptons, 0.5);

        // We don't care about events with no isolated jets
        if (goodjets.empty()) {
            MSG_DEBUG("No jets in event");
            vetoEvent;
        }

        Jets jc_final;
        Jets jb_final;

        // identification of bjets
        int n_btag = count(goodjets, hasBTag());

        for (const Jet& j : goodjets) {
            if (j.cTagged() && n_btag == 0) {
                jc_final.push_back(j);
            }
            if (j.bTagged()) {
                jb_final.push_back(j);
            }
        }
        // histogram filling

        if (goodjets.size() > 0) {
            if (jc_final.size() > 0) {
                _h_pt_cjet->fill(jc_final[0].pt());
                _h_Z_pt_cjet->fill(zcand.pt());
            }
            if (jb_final.size() > 0) {
                _h_pt_bjet->fill(jb_final[0].pt());
                _h_Z_pt_bjet->fill(zcand.pt());
            }
        }
    }

    /// Normalise histograms etc., after the run
    void finalize() {
        double norm = (sumOfWeights() != 0) ? crossSection() / picobarn / sumOfWeights() : 1.0;
        // account if we have electrons and muons in the sample
        if ((counter_ee->val() > 1.) && (counter_mm->val() > 1.)) {
            norm = norm / 2.;
        }

        divide(_h_pt_cjet, _h_pt_bjet, _h_R_jet_pt_cjet);
        divide(_h_Z_pt_cjet, _h_Z_pt_bjet, _h_R_Z_pt_cjet);

        scale(_h_Z_pt_cjet, norm);
        scale(_h_pt_cjet, norm);
    }

    ///@}

    /// @name Histograms

    Histo1DPtr _h_pt_cjet;
    Histo1DPtr _h_pt_bjet;
    Histo1DPtr _h_Z_pt_cjet, _h_Z_pt_bjet;

    Scatter2DPtr _h_R_jet_pt_cjet;
    Scatter2DPtr _h_R_Z_pt_cjet;

    CounterPtr counter_ee, counter_mm;
};

RIVET_DECLARE_PLUGIN(CMS_2017_I1634835);

}  // namespace Rivet