rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2017_I1499471

Measurements of the associated production of a Z boson and b jets in pp collisions at $\sqrt{s} = 8$ TeV
Experiment: CMS (LHC)
Inspire ID: 1499471
Status: VALIDATED
Authors:
  • Fabio Cossutti
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • Run MC generators with Z decaying leptonically into both electrons and muons at 8 TeV CoM energy. For jet-flavour inclusive samples, order of 20 million unweighted events can give a reasonable global comparison, but precision in the high jet multiplicity region/high jet pt may require substantially larger samples or statistical enhancement of high jet multiplicities.

Measurements of the associated production of a Z boson with at least one jet originating from a b quark in proton-proton collisions at $\sqrt{s} = 8$ TeV are presented. Differential cross sections are measured with data collected by the CMS experiment corresponding to an integrated luminosity of 19.8 $\text{fb}^{-1}$. Z bosons are reconstructed through their decays to electrons and muons. Cross sections are measured as a function of observables characterizing the kinematics of the b jet and the Z boson. Ratios of differential cross sections for the associated production with at least one b jet to the associated production with any jet are also presented. The production of a Z boson with at least two b jets is investigated, and differential cross sections are measured for the dijet system. Results are compared to theoretical predictions, testing two different flavour schemes for the choice of initial-state partons. Note of the authors: the Z boson selection and b-tagging algorithms used in the code are based on standard Rivet code and do not exactly coincide with those used in the paper. They anyway reproduce the paper results at the percent level or better for most of the bins.

Source code: CMS_2017_I1499471.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Tools/Cuts.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/VetoedFinalState.hh"
  5#include "Rivet/Projections/LeptonFinder.hh"
  6#include "Rivet/Projections/FastJets.hh"
  7#include "Rivet/Projections/DileptonFinder.hh"
  8
  9namespace Rivet {
 10
 11
 12  class CMS_2017_I1499471 : public Analysis {
 13  public:
 14
 15    /// Constructor
 16    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1499471);
 17
 18
 19    /// Book histograms and initialise projections before the run
 20    void init() {
 21
 22      DileptonFinder zeeFinder(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV &&
 23                               Cuts::abspid == PID::ELECTRON, Cuts::massIn(71.0*GeV, 111.0*GeV));
 24      declare(zeeFinder, "ZeeFinder");
 25
 26      DileptonFinder zmumuFinder(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV &&
 27                                 Cuts::abspid == PID::MUON, Cuts::massIn(71.0*GeV, 111.0*GeV));
 28      declare(zmumuFinder, "ZmumuFinder");
 29
 30      VisibleFinalState visfs;
 31      VetoedFinalState jetConstits(visfs);
 32      jetConstits.addVetoOnThisFinalState(zeeFinder);
 33      jetConstits.addVetoOnThisFinalState(zmumuFinder);
 34
 35      FastJets akt05Jets(jetConstits, JetAlg::ANTIKT, 0.5);
 36      declare(akt05Jets, "AntiKt05Jets");
 37
 38      //Histograms booking
 39
 40      book(_h_first_bjet_pt_b ,1,1,1);
 41      book(_h_first_bjet_abseta_b ,3,1,1);
 42      book(_h_Z_pt_b ,5,1,1);
 43      book(_h_HT_b ,7,1,1);
 44      book(_h_Dphi_Zb_b ,9,1,1);
 45
 46      book(_h_first_jet_pt_ratio ,2,1,1);
 47      book(_h_first_jet_abseta_ratio ,4,1,1);
 48      book(_h_Z_pt_ratio ,6,1,1);
 49      book(_h_HT_ratio ,8,1,1);
 50      book(_h_Dphi_Zj_ratio ,10,1,1);
 51
 52      book(_h_first_jet_pt, "first_jet_pt", refData(1,1,1) ); // (*_h_first_bjet_pt_b);
 53      book(_h_first_jet_abseta, "first_jet_abseta", refData(3,1,1) ); // (*_h_first_bjet_abseta_b);
 54      book(_h_Z_pt, "Z_pt", refData(5,1,1) ); // (*_h_Z_pt_b);
 55      book(_h_HT, "HT", refData(7,1,1) ); // (*_h_HT_b);
 56      book(_h_Dphi_Zj, "Dphi_Zj", refData(9,1,1) ); // (*_h_Dphi_Zb_b);
 57
 58      book(_h_first_bjet_pt_bb ,11,1,1);
 59      book(_h_second_bjet_pt_bb ,12,1,1);
 60      book(_h_Z_pt_bb ,13,1,1);
 61      book(_h_bb_mass_bb ,14,1,1);
 62      book(_h_Zbb_mass_bb ,15,1,1);
 63      book(_h_Dphi_bb ,16,1,1);
 64      book(_h_DR_bb ,17,1,1);
 65      book(_h_DR_Zbmin_bb ,18,1,1);
 66      book(_h_A_DR_Zb_bb ,19,1,1);
 67
 68      book(_h_bjet_multiplicity ,20,1,1);
 69
 70    }
 71
 72
 73    /// Perform the per-event analysis
 74    void analyze(const Event& event) {
 75
 76      const DileptonFinder& zeeFS = apply<DileptonFinder>(event, "ZeeFinder");
 77      const DileptonFinder& zmumuFS = apply<DileptonFinder>(event, "ZmumuFinder");
 78
 79      const Particles& zees = zeeFS.bosons();
 80      const Particles& zmumus = zmumuFS.bosons();
 81
 82      // We did not find exactly one Z. No good.
 83      if (zees.size() + zmumus.size() != 1) {
 84        MSG_DEBUG("Did not find exactly one good Z candidate");
 85        vetoEvent;
 86      }
 87
 88      //event identification depending on mass window
 89      bool ee_event=false;
 90      bool mm_event=false;
 91
 92      if (zees.size() == 1) { ee_event = true; }
 93      if (zmumus.size() == 1) { mm_event = true; }
 94      const Particles& theLeptons = zees.size() ? zeeFS.constituents() : zmumuFS.constituents();
 95
 96      // Cluster jets
 97      // NB. Veto has already been applied on leptons and photons used for dressing
 98      const FastJets& fj = apply<FastJets>(event, "AntiKt05Jets");
 99      const Jets& jets = fj.jetsByPt(Cuts::abseta < 2.4 && Cuts::pT > 30*GeV);
100
101      // Perform lepton-jet overlap and HT calculation
102      double Ht = 0;
103      Jets goodjets;
104      for (const Jet& j : jets) {
105        // Decide if this jet is "good", i.e. isolated from the leptons
106        /// @todo Nice use-case for any() and a C++11 lambda
107        bool overlap = false;
108        for (const Particle& l : theLeptons) {
109          if (Rivet::deltaR(j, l) < 0.5) {
110            overlap = true;
111            break;
112          }
113        }
114
115        // Fill HT and good-jets collection
116        if (overlap) continue;
117        goodjets.push_back(j);
118        Ht += j.pT();
119      }
120
121      // We don't care about events with no isolated jets
122      if (goodjets.empty()) {
123        MSG_DEBUG("No jets in event");
124        vetoEvent;
125      }
126
127      Jets jb_final;
128
129      //identification of bjets
130
131      for (const Jet& j : goodjets) {
132        if ( j.bTagged() ) { jb_final.push_back(j); }
133      }
134
135      //Event weight
136      const double w = 0.5;
137
138      //histogram filling
139
140      if ((ee_event || mm_event) && goodjets.size() > 0) {
141
142        FourMomentum j1(goodjets[0].momentum());
143
144        _h_first_jet_pt->fill(j1.pt(),w);
145        _h_first_jet_abseta->fill(fabs(j1.eta()),w);
146        if ( ee_event ) _h_Z_pt->fill(zees[0].pt(),w);
147        if ( mm_event ) _h_Z_pt->fill(zmumus[0].pt(),w);
148        _h_HT->fill(Ht,w);
149        if ( ee_event ) _h_Dphi_Zj->fill(deltaPhi(zees[0], j1),w);
150        if ( mm_event ) _h_Dphi_Zj->fill(deltaPhi(zmumus[0], j1),w);
151
152        if ( jb_final.size() > 0 ) {
153
154          FourMomentum b1(jb_final[0].momentum());
155
156          _h_bjet_multiplicity->fill(1.,w);
157
158          _h_first_bjet_pt_b->fill(b1.pt(),w);
159          _h_first_bjet_abseta_b->fill(fabs(b1.eta()),w);
160          if ( ee_event ) _h_Z_pt_b->fill(zees[0].pt(),w);
161          if ( mm_event ) _h_Z_pt_b->fill(zmumus[0].pt(),w);
162          _h_HT_b->fill(Ht,w);
163          if ( ee_event ) _h_Dphi_Zb_b->fill(deltaPhi(zees[0], b1.phi()),w);
164          if ( mm_event ) _h_Dphi_Zb_b->fill(deltaPhi(zmumus[0], b1.phi()),w);
165
166          if ( jb_final.size() > 1 ) {
167
168            FourMomentum b2(jb_final[1].momentum());
169
170            _h_bjet_multiplicity->fill(2.,w);
171
172            _h_first_bjet_pt_bb->fill(b1.pt(),w);
173            _h_second_bjet_pt_bb->fill(b2.pt(),w);
174            if ( ee_event ) _h_Z_pt_bb->fill(zees[0].pt(),w);
175            if ( mm_event ) _h_Z_pt_bb->fill(zmumus[0].pt(),w);
176
177            FourMomentum bb = add(b1,b2);
178            FourMomentum Zbb;
179            if (ee_event) Zbb = add(zees[0],bb);
180            if (mm_event) Zbb = add(zmumus[0],bb);
181
182            _h_bb_mass_bb->fill(bb.mass(),w);
183            _h_Zbb_mass_bb->fill(Zbb.mass(),w);
184
185            _h_Dphi_bb->fill(deltaPhi(b1,b2),w);
186            if (deltaR(b1,b2)>0.5) {
187              _h_DR_bb->fill(deltaR(b1,b2),w);
188            }
189
190            double DR_Z_b1(0.), DR_Z_b2(0.);
191            if ( ee_event ) {
192              DR_Z_b1 = deltaR(zees[0],b1);
193              DR_Z_b2 = deltaR(zees[0],b2);
194            }
195            if ( mm_event ) {
196              DR_Z_b1 = deltaR(zmumus[0],b1);
197              DR_Z_b2 = deltaR(zmumus[0],b2);
198            }
199
200            double DR_Zb_min = DR_Z_b1;
201            double DR_Zb_max = DR_Z_b2;
202            if ( DR_Zb_min > DR_Zb_max ) {
203              DR_Zb_min = DR_Z_b2;
204              DR_Zb_max = DR_Z_b1;
205            }
206            double A_Zbb = (DR_Zb_max - DR_Zb_min)/(DR_Zb_max + DR_Zb_min);
207
208            _h_DR_Zbmin_bb->fill(DR_Zb_min,w);
209            _h_A_DR_Zb_bb->fill(A_Zbb,w);
210
211          }
212
213        }
214
215      }
216
217    }
218
219
220    /// Normalise histograms etc., after the run
221    void finalize() {
222
223      const double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/picobarn/sumOfWeights() : 1.0;
224
225      MSG_INFO("Cross section = " << std::setfill(' ') << std::setw(14) << std::fixed << std::setprecision(3) << crossSection()/picobarn << " pb");
226      MSG_INFO("# Events      = " << std::setfill(' ') << std::setw(14) << std::fixed << std::setprecision(3) << numEvents() );
227      MSG_INFO("SumW          = " << std::setfill(' ') << std::setw(14) << std::fixed << std::setprecision(3) << sumOfWeights());
228      MSG_INFO("Norm factor   = " << std::setfill(' ') << std::setw(14) << std::fixed << std::setprecision(6) << norm);
229
230      scale( _h_first_bjet_pt_b, 100. );
231      scale( _h_first_bjet_abseta_b, 100. );
232      scale( _h_Z_pt_b, 100. );
233      scale( _h_HT_b, 100. );
234      scale( _h_Dphi_Zb_b, 100. );
235
236      divide( _h_first_bjet_pt_b , _h_first_jet_pt , _h_first_jet_pt_ratio );
237      divide( _h_first_bjet_abseta_b , _h_first_jet_abseta , _h_first_jet_abseta_ratio );
238      divide( _h_Z_pt_b , _h_Z_pt , _h_Z_pt_ratio );
239      divide( _h_HT_b , _h_HT , _h_HT_ratio );
240      divide( _h_Dphi_Zb_b , _h_Dphi_Zj , _h_Dphi_Zj_ratio );
241
242      scale( _h_first_bjet_pt_b, norm/100. );
243      scale( _h_first_bjet_abseta_b, norm/100. );
244      scale( _h_Z_pt_b, norm/100. );
245      scale( _h_HT_b, norm/100. );
246      scale( _h_Dphi_Zb_b, norm/100. );
247
248      scale( _h_first_bjet_pt_bb, norm);
249      scale( _h_second_bjet_pt_bb, norm);
250      scale( _h_Z_pt_bb, norm);
251      scale( _h_bb_mass_bb, norm);
252      scale( _h_Zbb_mass_bb, norm);
253      scale( _h_Dphi_bb, norm);
254      scale( _h_DR_bb, norm);
255      scale( _h_DR_Zbmin_bb, norm);
256      scale( _h_A_DR_Zb_bb, norm);
257
258      scale( _h_bjet_multiplicity, norm );
259
260    }
261
262
263  private:
264
265    /// @name Histograms
266
267    Histo1DPtr     _h_first_jet_pt, _h_first_bjet_pt_b;
268    Histo1DPtr     _h_first_jet_abseta, _h_first_bjet_abseta_b;
269    Histo1DPtr     _h_Z_pt, _h_Z_pt_b;
270    Histo1DPtr     _h_HT, _h_HT_b;
271    Histo1DPtr     _h_Dphi_Zj, _h_Dphi_Zb_b;
272
273    Estimate1DPtr     _h_first_jet_pt_ratio;
274    Estimate1DPtr     _h_first_jet_abseta_ratio;
275    Estimate1DPtr     _h_Z_pt_ratio;
276    Estimate1DPtr     _h_HT_ratio;
277    Estimate1DPtr     _h_Dphi_Zj_ratio;
278
279    Histo1DPtr     _h_first_bjet_pt_bb, _h_second_bjet_pt_bb;
280    Histo1DPtr     _h_Z_pt_bb;
281    Histo1DPtr     _h_bb_mass_bb, _h_Zbb_mass_bb;
282    Histo1DPtr     _h_Dphi_bb, _h_DR_bb, _h_DR_Zbmin_bb, _h_A_DR_Zb_bb;
283
284    Histo1DPtr     _h_bjet_multiplicity;
285
286  };
287
288
289  RIVET_DECLARE_PLUGIN(CMS_2017_I1499471);
290
291}