rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1469071

Measurement of the $WZ$ production cross section at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1469071
Status: VALIDATED
Authors:
  • Elena Yatsenko
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • pp -> WZ + X, diboson decays to electrons and muons

The production of $W^\pm Z$ events in proton--proton collisions at a centre-of-mass energy of 13 TeV is measured with the ATLAS detector at the LHC. The collected data correspond to an integrated luminosity of 3.2 fb${}^{-1}$ . The $W^\pm Z$ candidates are reconstructed using leptonic decays of the gauge bosons into electrons or muons. The measured inclusive cross section in the detector fiducial region for leptonic decay modes is $\sigma(W^\pm Z\to\ell^\prime\nu\ell\ell \text{fid.}) = 63.2\pm 3.2$(stat.)$\pm 2.6$(sys.)$\pm 1.5$(lumi.) fb. In comparison, the next-to-leading-order Standard Model prediction is 53.4-2.8+3.6 fb. The extrapolation of the measurement from the fiducial to the total phase space yields $\sigma(W^\pm Z\text{tot.})=50.6\pm 2.6$(stat.)$\pm 2.0$(sys.)$\pm 0.9$(th.)$\pm 1.2$(lumi.) pb, in agreement with a recent next-to-next-to-leading-order calculation of 48.2-1.0+1.1 pb. The cross section as a function of jet multiplicity is also measured, together with the charge-dependent $W^+Z$ and $W^-Z$ cross sections and their ratio. Uses SM neutrino-lepton flavour matching and a resonant shape algorithm assuming the Standard Model, to match the MC-based correction to the fiducial region applied in the paper. This routine is therefore only valid under the assumption of the Standard Model and cannot be used for BSM reinterpretation.

Source code: ATLAS_2016_I1469071.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/IdentifiedFinalState.hh"
  6#include "Rivet/Projections/PromptFinalState.hh"
  7#include "Rivet/Projections/LeptonFinder.hh"
  8#include "Rivet/Projections/VetoedFinalState.hh"
  9
 10namespace Rivet {
 11
 12
 13  /// @brief Measurement of the WZ production cross section at 13 TeV
 14  class ATLAS_2016_I1469071 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1469071);
 19
 20    /// @name Analysis methods
 21    /// @{
 22
 23    /// Book histograms and initialise projections before the run
 24    void init() {
 25
 26      // Lepton cuts
 27      Cut FS_Zlept = Cuts::abseta < 2.5 && Cuts::pT > 15*GeV;
 28
 29      FinalState fs;
 30      Cut fs_z = Cuts::abseta < 2.5 && Cuts::pT > 15*GeV;
 31      Cut fs_j = Cuts::abseta < 4.5 && Cuts::pT > 25*GeV;
 32
 33      // Get photons to dress leptons
 34      PromptFinalState photons(Cuts::abspid == PID::PHOTON);
 35
 36      // Electrons and muons in Fiducial PS
 37      PromptFinalState leptons(FinalState(fs_z && (Cuts::abspid == PID::ELECTRON || Cuts::abspid == PID::MUON)));
 38      leptons.acceptTauDecays(false);
 39      LeptonFinder dressedleptons(leptons, photons, 0.1, FS_Zlept);
 40      declare(dressedleptons, "LeptonFinder");
 41
 42      // Electrons and muons in Total PS
 43      PromptFinalState leptons_total(Cuts::abspid == PID::ELECTRON || Cuts::abspid == PID::MUON);
 44      leptons_total.acceptTauDecays(false);
 45      LeptonFinder dressedleptonsTotal(leptons_total, photons, 0.1);
 46      declare(dressedleptonsTotal, "LeptonFinderTotal");
 47
 48      // Promot neutrinos (yikes!)
 49      IdentifiedFinalState nu_id;
 50      nu_id.acceptNeutrinos();
 51      PromptFinalState neutrinos(nu_id);
 52      neutrinos.acceptTauDecays(false);
 53      declare(neutrinos, "Neutrinos");
 54      MSG_WARNING("\033[91;1mLIMITED VALIDITY - check info file for details!\033[m");
 55
 56      // Jets
 57      VetoedFinalState veto;
 58      veto.addVetoOnThisFinalState(dressedleptons);
 59      FastJets jets(veto, JetAlg::ANTIKT, 0.4);
 60      declare(jets, "Jets");
 61
 62      // Book histograms
 63      book(_h_eee      , 1, 1, 1);
 64      book(_h_mee      , 1, 1, 2);
 65      book(_h_emm      , 1, 1, 3);
 66      book(_h_mmm      , 1, 1, 4);
 67      book(_h_fid      , 1, 1, 5);
 68      book(_h_eee_Plus , 2, 1, 1);
 69      book(_h_mee_Plus , 2, 1, 2);
 70      book(_h_emm_Plus , 2, 1, 3);
 71      book(_h_mmm_Plus , 2, 1, 4);
 72      book(_h_fid_Plus , 2, 1, 5);
 73      book(_h_eee_Minus, 3, 1, 1);
 74      book(_h_mee_Minus, 3, 1, 2);
 75      book(_h_emm_Minus, 3, 1, 3);
 76      book(_h_mmm_Minus, 3, 1, 4);
 77      book(_h_fid_Minus, 3, 1, 5);
 78      book(_h_total    , 6, 1, 1);
 79      book(_h_Njets    , 8, 1, 1);
 80
 81    }
 82
 83
 84    void analyze(const Event& event) {
 85
 86      const DressedLeptons& dressedleptons = apply<LeptonFinder>(event, "LeptonFinder").dressedLeptons();
 87      const DressedLeptons& dressedleptonsTotal = apply<LeptonFinder>(event, "LeptonFinderTotal").dressedLeptons();
 88      const Particles& neutrinos = apply<PromptFinalState>(event, "Neutrinos").particlesByPt();
 89      Jets jets = apply<JetFinder>(event, "Jets").jetsByPt( (Cuts::abseta < 4.5) && (Cuts::pT > 25*GeV) );
 90
 91      if (dressedleptonsTotal.size() < 3 || neutrinos.size() < 1) vetoEvent;
 92
 93      //---Total PS: assign leptons to W and Z bosons using Resonant shape algorithm
 94      // NB: This resonant shape algorithm assumes the Standard Model and can therefore
 95      //     NOT be used for any kind of reinterpretation in terms of new-physics models..
 96
 97      // Try Z pair of leptons 01
 98      double MassZ01 = 0, MassW2 = 0;
 99      if ( (dressedleptonsTotal[0].pid() ==-(dressedleptonsTotal[1].pid())) && (dressedleptonsTotal[2].abspid()==neutrinos[0].abspid()-1)){
100        MassZ01 = (dressedleptonsTotal[0].momentum()+dressedleptonsTotal[1].momentum()).mass();
101        MassW2 = (dressedleptonsTotal[2].momentum()+neutrinos[0].momentum()).mass();
102      }
103      // Try Z pair of leptons 02
104      double MassZ02 = 0, MassW1 = 0;
105      if ( (dressedleptonsTotal[0].pid()==-(dressedleptonsTotal[2].pid())) && (dressedleptonsTotal[1].abspid()==neutrinos[0].abspid()-1)){
106        MassZ02 = (dressedleptonsTotal[0].momentum()+dressedleptonsTotal[2].momentum()).mass();
107        MassW1 = (dressedleptonsTotal[1].momentum()+neutrinos[0].momentum()).mass();
108      }
109      // Try Z pair of leptons 12
110      double MassZ12 = 0, MassW0 = 0;
111      if ( (dressedleptonsTotal[1].pid()==-(dressedleptonsTotal[2].pid())) && (dressedleptonsTotal[0].abspid()==neutrinos[0].abspid()-1)){
112        MassZ12 = (dressedleptonsTotal[1].momentum()+dressedleptonsTotal[2].momentum()).mass();
113        MassW0 = (dressedleptonsTotal[0].momentum()+neutrinos[0].momentum()).mass();
114      }
115      double WeightZ1 = 1/(pow(MassZ01*MassZ01 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
116      double WeightW1 = 1/(pow(MassW2*MassW2 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
117      double WeightTotal1 = WeightZ1*WeightW1;
118      double M1 = -1*WeightTotal1;
119
120      double WeightZ2 = 1/(pow(MassZ02*MassZ02- MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
121      double WeightW2 = 1/(pow(MassW1*MassW1- MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
122      double WeightTotal2 = WeightZ2*WeightW2;
123      double M2 = -1*WeightTotal2;
124
125      double WeightZ3 = 1/(pow(MassZ12*MassZ12 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
126      double WeightW3 = 1/(pow(MassW0*MassW0 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
127      double WeightTotal3 = WeightZ3*WeightW3;
128      double M3 = -1*WeightTotal3;
129
130      int i = -1, j = -1, k = -1;
131      if ((M1 < M2 && M1 < M3) || (MassZ01 != 0 && MassW2 != 0 && MassZ02 == 0 && MassZ12 == 0)) {
132        i = 0; j = 1; k = 2;
133      }
134      if ((M2 < M1 && M2 < M3) || (MassZ02 != 0 && MassW1 != 0 && MassZ01 == 0 && MassZ12 == 0)) {
135        i = 0; j = 2; k = 1;
136      }
137      if ((M3 < M1 && M3 < M2) || (MassZ12 != 0 && MassW0 != 0 && MassZ01 == 0 && MassZ02 == 0)) {
138        i = 1; j = 2; k = 0;
139      }
140      ///
141      if (i < 0 || j < 0 || k < 0) vetoEvent;
142
143      FourMomentum ZbosonTotal = dressedleptonsTotal[i].momentum()+dressedleptonsTotal[j].momentum();
144      if ( ZbosonTotal.mass() >= 66*GeV && ZbosonTotal.mass() <= 116*GeV )  _h_total->fill(13000);
145
146      //---end Total PS
147
148
149      //---Fiducial PS: assign leptons to W and Z bosons using Resonant shape algorithm
150      if (dressedleptons.size() < 3)  vetoEvent;
151
152      int EventType = -1;
153      int Nel = 0, Nmu = 0;
154
155      for (const DressedLepton& l : dressedleptons) {
156        if (l.abspid() == 11)  ++Nel;
157        if (l.abspid() == 13)  ++Nmu;
158      }
159
160      if ( (Nel == 3)  && (Nmu==0) )  EventType = 3;
161      if ( (Nel == 2)  && (Nmu==1) )  EventType = 2;
162      if ( (Nel == 1)  && (Nmu==2) )  EventType = 1;
163      if ( (Nel == 0)  && (Nmu==3) )  EventType = 0;
164
165      int EventCharge = -dressedleptons[0].charge() * dressedleptons[1].charge() * dressedleptons[2].charge();
166
167      MassZ01 = 0; MassZ02 = 0; MassZ12 = 0;
168      MassW0 = 0;  MassW1 = 0;  MassW2 = 0;
169
170      // try Z pair of leptons 01
171      if (dressedleptons[0].pid() == -dressedleptons[1].pid()) {
172        MassZ01 = (dressedleptons[0].momentum() + dressedleptons[1].momentum()).mass();
173        MassW2 = (dressedleptons[2].momentum() + neutrinos[0].momentum()).mass();
174      }
175      // try Z pair of leptons 02
176      if (dressedleptons[0].pid() == -dressedleptons[2].pid()) {
177        MassZ02 = (dressedleptons[0].momentum() + dressedleptons[2].momentum()).mass();
178        MassW1 = (dressedleptons[1].momentum() + neutrinos[0].momentum()).mass();
179      }
180      // try Z pair of leptons 12
181      if (dressedleptons[1].pid() == -dressedleptons[2].pid()) {
182        MassZ12 = (dressedleptons[1].momentum() + dressedleptons[2].momentum()).mass();
183        MassW0 = (dressedleptons[0].momentum() + neutrinos[0].momentum()).mass();
184      }
185      WeightZ1 = 1/(pow(MassZ01*MassZ01 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
186      WeightW1 = 1/(pow(MassW2*MassW2 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
187      WeightTotal1 = WeightZ1*WeightW1;
188      M1 = -1*WeightTotal1;
189
190      WeightZ2 = 1/(pow(MassZ02*MassZ02- MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
191      WeightW2 = 1/(pow(MassW1*MassW1- MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
192      WeightTotal2 = WeightZ2*WeightW2;
193      M2 = -1*WeightTotal2;
194
195      WeightZ3 = 1/(pow(MassZ12*MassZ12 - MZ_PDG*MZ_PDG,2) + pow(MZ_PDG*GammaZ_PDG,2));
196      WeightW3 = 1/(pow(MassW0*MassW0 - MW_PDG*MW_PDG,2) + pow(MW_PDG*GammaW_PDG,2));
197      WeightTotal3 = WeightZ3*WeightW3;
198      M3 = -1*WeightTotal3;
199
200      if( (M1 < M2 && M1 < M3) || (MassZ01 != 0 && MassW2 != 0 && MassZ02 == 0 && MassZ12 == 0) ) {
201        i = 0; j = 1; k = 2;
202      }
203      if((M2 < M1 && M2 < M3) || (MassZ02 != 0 && MassW1 != 0 && MassZ01 == 0 && MassZ12 == 0) ) {
204        i = 0; j = 2; k = 1;
205      }
206      if((M3 < M1 && M3 < M2) || (MassZ12 != 0 && MassW0 != 0 && MassZ01 == 0 && MassZ02 == 0) ) {
207        i = 1; j = 2; k = 0;
208      }
209
210      FourMomentum Zlepton1 = dressedleptons[i].momentum();
211      FourMomentum Zlepton2 = dressedleptons[j].momentum();
212      FourMomentum Wlepton  = dressedleptons[k].momentum();
213      FourMomentum Zboson   = dressedleptons[i].momentum()+dressedleptons[j].momentum();
214
215      double Wboson_mT = sqrt( 2 * Wlepton.pT() * neutrinos[0].pt() * (1 - cos(deltaPhi(Wlepton, neutrinos[0]))) );
216
217      if (fabs(Zboson.mass()/GeV - MZ_PDG) >= 10.) vetoEvent;
218      if (Wboson_mT <= 30*GeV)                     vetoEvent;
219      if (Wlepton.pT() <= 20*GeV)                  vetoEvent;
220      if (deltaR(Zlepton1, Zlepton2) < 0.2)        vetoEvent;
221      if (deltaR(Zlepton1, Wlepton)  < 0.3)        vetoEvent;
222      if (deltaR(Zlepton2, Wlepton)  < 0.3)        vetoEvent;
223
224      if (EventType == 3)  _h_eee->fill(13000.);
225      if (EventType == 2)  _h_mee->fill(13000.);
226      if (EventType == 1)  _h_emm->fill(13000.);
227      if (EventType == 0)  _h_mmm->fill(13000.);
228      _h_fid->fill(13000);
229
230      if (EventCharge == 1) {
231        if (EventType == 3)  _h_eee_Plus->fill(13000.);
232        if (EventType == 2)  _h_mee_Plus->fill(13000.);
233        if (EventType == 1)  _h_emm_Plus->fill(13000.);
234        if (EventType == 0)  _h_mmm_Plus->fill(13000.);
235        _h_fid_Plus->fill(13000);
236      } else {
237        if (EventType == 3)  _h_eee_Minus->fill(13000.);
238        if (EventType == 2)  _h_mee_Minus->fill(13000.);
239        if (EventType == 1)  _h_emm_Minus->fill(13000.);
240        if (EventType == 0)  _h_mmm_Minus->fill(13000.);
241        _h_fid_Minus->fill(13000);
242      }
243
244      if (jets.size() < 4)  _h_Njets->fill(jets.size());
245      else  _h_Njets->fill(4);
246
247    }
248
249
250    void finalize() {
251
252      // Print summary info
253      const double xs_pb(crossSection() / picobarn);
254      const double xs_fb(crossSection() / femtobarn);
255      const double sumw(sumOfWeights());
256      const double sf_pb(xs_pb / sumw);
257      const double sf_fb(xs_fb / sumw);
258
259      const float totalBR= 4*0.1086*0.033658; // W and Z leptonic branching fractions
260
261      scale(_h_fid,       sf_fb/4.);
262      scale(_h_eee,       sf_fb);
263      scale(_h_mee,       sf_fb);
264      scale(_h_emm,       sf_fb);
265      scale(_h_mmm,       sf_fb);
266      scale(_h_fid_Plus,  sf_fb/4.);
267      scale(_h_eee_Plus,  sf_fb);
268      scale(_h_mee_Plus,  sf_fb);
269      scale(_h_emm_Plus,  sf_fb);
270      scale(_h_mmm_Plus,  sf_fb);
271      scale(_h_fid_Minus, sf_fb/4.);
272      scale(_h_eee_Minus, sf_fb);
273      scale(_h_mee_Minus, sf_fb);
274      scale(_h_emm_Minus, sf_fb);
275      scale(_h_mmm_Minus, sf_fb);
276      scale(_h_Njets, sf_fb/4.);
277      scale(_h_total, sf_pb/totalBR);
278
279    }
280
281    /// @}
282
283
284  private:
285
286
287    /// @name Histograms
288    /// @{
289    Histo1DPtr _h_eee;
290    Histo1DPtr _h_mee;
291    Histo1DPtr _h_emm;
292    Histo1DPtr _h_mmm;
293    Histo1DPtr _h_fid;
294    Histo1DPtr _h_eee_Plus;
295    Histo1DPtr _h_mee_Plus;
296    Histo1DPtr _h_emm_Plus;
297    Histo1DPtr _h_mmm_Plus;
298    Histo1DPtr _h_fid_Plus;
299    Histo1DPtr _h_eee_Minus;
300    Histo1DPtr _h_mee_Minus;
301    Histo1DPtr _h_emm_Minus;
302    Histo1DPtr _h_mmm_Minus;
303    Histo1DPtr _h_fid_Minus;
304    Histo1DPtr _h_total;
305    Histo1DPtr _h_Njets;
306
307    /// @}
308
309    double MZ_PDG = 91.1876;
310    double MW_PDG = 83.385;
311    double GammaZ_PDG = 2.4952;
312    double GammaW_PDG = 2.085;
313
314  };
315
316  RIVET_DECLARE_PLUGIN(ATLAS_2016_I1469071);
317}