rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1448301

$Z \gamma (\gamma)$ cross sections at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1448301
Status: VALIDATED
Authors:
  • Hulin Wang
  • Evgeny Soldatov
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • p + p -> e+ e- gamma (gamma) at 8 TeV

The production of $Z$ bosons with one or two isolated high-energy photons is studied using pp collisions at $\sqrt{s}=8$ TeV. The analyses use a data sample with an integrated luminosity of 20.3 fb${}^{-1}$ collected by the ATLAS detector during the 2012 LHC data taking. The $Z\gamma$ and $Z\gamma\gamma$ production cross sections are measured with leptonic ($e^+e^-$, $\mu^+\mu^-$, $\nu\bar{\nu}$) decays of the Z boson, in extended fiducial regions defined in terms of the lepton and photon acceptance. They are then compared to cross-section predictions from the Standard Model, where the sources of the photons are radiation off initial-state quarks and radiative Z-boson decay to charged leptons, and from fragmentation of final-state quarks and gluons into photons. The yields of events with photon transverse energy $E_\text{T}>250$ GeV from $\ell^+\ell^-\gamma$ events and with $E_\text{T}>400$ GeV from $\nu\bar{\nu}\gamma$ events are used to search for anomalous triple gauge-boson couplings $ZZ\gamma$ and $Z\gamma\gamma$. The yields of events with diphoton invariant mass $m_{\gamma\gamma}>200$ GeV from $\ell^+\ell^-\gamma\gamma$ events and with $m_{\gamma\gamma}>300$ GeV from $\nu\bar{\nu}\gamma\gamma$ events are used to search for anomalous quartic gauge-boson couplings $ZZ\gamma\gamma$ and $Z\gamma\gamma\gamma$. No deviations from Standard Model predictions are observed and limits are placed on parameters used to describe anomalous triple and quartic gauge-boson couplings. The lepton channel can be specified by using the relevant LMODE (EL, MU, NU, LL). The default plugin will fill all of them with any events passing the cuts. NU will do only the MET channels, EL & MU will do only the electron or muon and will populate the combined plots assuming charged lepton universality. LL with do on the electron and muon but will populate the combined plots with the average of the two.

Source code: ATLAS_2016_I1448301.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/VetoedFinalState.hh"
  6#include "Rivet/Projections/IdentifiedFinalState.hh"
  7#include "Rivet/Projections/LeptonFinder.hh"
  8#include "Rivet/Projections/PromptFinalState.hh"
  9#include "Rivet/Projections/VisibleFinalState.hh"
 10
 11namespace Rivet {
 12
 13
 14  /// @brief Z/gamma cross section measurement at 8 TeV
 15  class ATLAS_2016_I1448301 : public Analysis {
 16  public:
 17
 18    /// Constructor
 19    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1448301);
 20
 21
 22    /// @name Analysis methods
 23    /// @{
 24
 25    /// Book histograms and initialise projections before the run
 26    void init() {
 27
 28      // Get options from the new option system
 29      // Default tries to fill everything
 30      // NU fills only the MET plots
 31      // EL fills only the Electron plots and combined plots assuming lepton univerality
 32      // MU fills only the Muon plots and combined plots assuming lepton univerality
 33      // LL fills electron and Muon plots and combined plots from correct average
 34      _mode = 0;
 35      if ( getOption("LMODE") == "NU" ) _mode = 1;
 36      if ( getOption("LMODE") == "EL" ) _mode = 2;
 37      if ( getOption("LMODE") == "MU" ) _mode = 3;
 38      if ( getOption("LMODE") == "LL" ) _mode = 4;
 39
 40      // Prompt photons
 41      const Cut photoncut = Cuts::abspid == PID::PHOTON && Cuts::pT > 15*GeV && Cuts::abseta < 2.37;
 42      PromptFinalState photon_fs(photoncut);
 43      declare(photon_fs, "Photons");
 44
 45      // Prompt leptons
 46      const PromptFinalState bareelectron_fs = Cuts::abspid == PID::ELECTRON;
 47      const PromptFinalState baremuon_fs = Cuts::abspid == PID::MUON;
 48
 49      // Dressed leptons
 50      const IdentifiedFinalState allphoton_fs(PID::PHOTON); // photons used for lepton dressing
 51      const Cut leptoncut = Cuts::pT > 25*GeV && Cuts::abseta < 2.47;
 52      const LeptonFinder dressedelectron_fs(bareelectron_fs, allphoton_fs, 0.1, leptoncut);
 53      const LeptonFinder dressedmuon_fs(baremuon_fs, allphoton_fs, 0.1, leptoncut);
 54
 55      declare(dressedelectron_fs, "Electrons");
 56      declare(dressedmuon_fs, "Muons");
 57
 58      // MET (prompt neutrinos)
 59      VetoedFinalState ivfs;
 60      ivfs.addVetoOnThisFinalState(VisibleFinalState());
 61      declare(PromptFinalState(ivfs), "MET");
 62
 63      // Jets
 64      VetoedFinalState jet_fs;
 65      jet_fs.vetoNeutrinos();
 66      jet_fs.addVetoPairId(PID::MUON);
 67      const FastJets fastjets(jet_fs, JetAlg::ANTIKT, 0.4);
 68      declare(fastjets, "Jets");
 69
 70      // Histograms
 71
 72      // MET
 73      if (_mode == 0 || _mode == 1){
 74        book(_h["vvg"],     2, 1, 1);
 75        book(_h["vvgg"],    4, 1, 1);
 76        book(_h["pT"],      7, 1, 1);
 77        book(_h["pT_0jet"], 8, 1, 1);
 78      }
 79
 80      // always book e and mu in charged lepton modes; there are sometimes 4 leptons.
 81      if (_mode != 1){
 82	// electron
 83        book(_h["eeg"],  1, 1, 1);
 84        book(_h["eegg"], 3, 1, 1);
 85        // muon
 86	book(_h["mmg"],  1, 1, 2);
 87	book(_h["mmgg"], 3, 1, 2);
 88
 89        // combined
 90        book(_h["llgg"], 3, 1, 3);
 91        book(_h["llg"], 1, 1, 3);
 92        book(_h["pT"], 5, 1, 1);
 93        book(_h["pT_0jet"], 6, 1, 1);
 94        book(_h["M"], 9, 1, 1);
 95        book(_h["M_0jet"], 10, 1, 1);
 96        book(_h["Njets"], 11, 1, 1);
 97      }
 98    }
 99
100
101    /// Perform the per-event analysis
102    void analyze(const Event& event) {
103
104      // Get objects
105      DressedLeptons electrons = apply<LeptonFinder>(event, "Electrons").dressedLeptons();
106      DressedLeptons muons = apply<LeptonFinder>(event, "Muons").dressedLeptons();
107
108      const Particles& photons = apply<PromptFinalState>(event, "Photons").particlesByPt();
109      const Jets jets = apply<FastJets>(event, "Jets").jetsByPt();
110
111      const FinalState& metfs = apply<PromptFinalState>(event, "MET");
112      Vector3 met_vec;
113      for (const Particle& p : metfs.particles()) met_vec += p.mom().perpVec();
114
115      if (met_vec.mod() >= 100*GeV && !photons.empty() && _mode < 2){
116
117	if (photons.size() > 1) { // nu nu y y
118
119	  bool yy_veto = false;
120	  yy_veto |= photons[0].pT() < 22*GeV;
121	  yy_veto |= photons[1].pT() < 22*GeV;
122	  yy_veto |= met_vec.mod() < 110*GeV;
123	  const double yyPhi = (photons[0].momentum() + photons[1].momentum()).phi();
124	  yy_veto |= fabs(yyPhi - met_vec.phi()) < 2.62 || fabs(yyPhi - met_vec.phi()) > 3.66;
125	  yy_veto |= deltaR(photons[0], photons[1]) < 0.4;
126
127	  // Photon isolation calculated by jets, count jets
128	  Jet ph0_jet, ph1_jet;
129	  double min_dR_ph0_jet = 999., min_dR_ph1_jet = 999.;
130	  size_t njets = 0;
131	  for (const Jet& j : jets) {
132	    if (j.pT() > 30*GeV && j.abseta() < 4.5) {
133	      if (deltaR(j, photons[0]) > 0.3 && deltaR(j, photons[1]) > 0.3)  ++njets;
134	    }
135	    if (deltaR(j, photons[0]) < min_dR_ph0_jet) {
136	      min_dR_ph0_jet = deltaR(j, photons[0]);
137	      ph0_jet = j;
138	    }
139	    if (deltaR(j, photons[1]) < min_dR_ph1_jet) {
140	      min_dR_ph1_jet = deltaR(j, photons[1]);
141	      ph1_jet = j;
142	    }
143	  }
144	  double photon0iso = 0., photon1iso = 0.;
145	  if (min_dR_ph0_jet < 0.4)  photon0iso = ph0_jet.pT() - photons[0].pT();
146	  if (min_dR_ph1_jet < 0.4)  photon1iso = ph1_jet.pT() - photons[1].pT();
147	  yy_veto |= photon0iso/photons[0].pT() > 0.5;
148	  yy_veto |= photon1iso/photons[1].pT() > 0.5;
149
150	  if (!yy_veto) {
151	    _h["vvgg"]->fill(0.5);
152	    if (!njets)  _h["vvgg"]->fill(1.5);
153	  }
154	} // end of nu nu y y section
155
156
157	if ((photons[0].pT() >= 130*GeV)  &&
158	    (fabs(fabs(deltaPhi(photons[0], met_vec)) - 3.14) <= 1.57)) {
159
160	  // Photon isolation calculated by jets, count jets
161	  Jet ph_jet;
162	  double min_dR_ph_jet = 999.;
163	  size_t njets = 0;
164	  for (const Jet& j : jets) {
165	    if (j.pT() > 30*GeV && j.abseta() < 4.5) {
166	      if (deltaR(j, photons[0]) > 0.3)  ++njets;
167	    }
168	    if (deltaR(j, photons[0]) < min_dR_ph_jet) {
169	      min_dR_ph_jet = deltaR(j, photons[0]);
170	      ph_jet = j;
171	    }
172	  }
173	  double photoniso = 0;
174	  if (min_dR_ph_jet < 0.4)  photoniso = ph_jet.pT() - photons[0].pT();
175	  if (photoniso/photons[0].pT() > 0.5)  vetoEvent;
176
177	  const double pTgamma = photons[0].pT()/GeV;
178	  _h["pT"]->fill(pTgamma);
179	  _h["vvg"]->fill(0.5);
180	  if (!njets) {
181	    _h["vvg"]->fill(1.5);
182	    _h["pT_0jet"]->fill(pTgamma);
183	  }
184
185	}
186      }  // end of nu nu y (y) section
187
188      // Dilepton candidate
189      bool el = false;
190      if ( (_mode != 1) &&
191	   (( electrons.size() >= 2 && _mode != 3 ) ||
192	    ( muons.size()     >= 2 && _mode != 2 ) )) {
193
194	DressedLeptons lep_p, lep_m;
195
196	// Sort the dressed leptons by pt
197	if (electrons.size() >= 2) {
198	  el = true;
199	  sortByPt(electrons);
200	  for (const DressedLepton& lep : electrons) {
201	    if (lep.charge() > 0.)  lep_p.push_back(lep);
202	    if (lep.charge() < 0.)  lep_m.push_back(lep);
203	  }
204	} else {
205	  sortByPt(muons);
206	  for (const DressedLepton& lep : muons) {
207	    if (lep.charge() > 0.)  lep_p.push_back(lep);
208	    if (lep.charge() < 0.)  lep_m.push_back(lep);
209	  }
210	}
211
212
213	if (!lep_p.empty() && !lep_m.empty() &&
214	    (lep_p[0].abspid() == lep_m[0].abspid()) &&
215	    ((lep_p[0].momentum() + lep_m[0].momentum()).mass() >= 40*GeV)){
216
217	  // Photon lepton overlap removal
218	  if (photons.empty())  vetoEvent;
219
220	  if (photons.size() > 1) {
221
222	    bool veto = false;
223	    veto |= deltaR(photons[0], lep_p[0]) < 0.4;
224	    veto |= deltaR(photons[0], lep_m[0]) < 0.4;
225	    veto |= deltaR(photons[1], lep_p[0]) < 0.4;
226	    veto |= deltaR(photons[1], lep_m[0]) < 0.4;
227	    veto |= deltaR(photons[0], photons[1]) < 0.4;
228
229	    Jet ph0_jet, ph1_jet;
230	    double min_dR_ph0_jet = 999., min_dR_ph1_jet=999.;
231	    int njets = 0;
232	    for (const Jet& j : jets){
233	      if (j.pT() > 30*GeV && j.abseta() < 4.5) {
234		if (deltaR(j, lep_p[0]) > 0.3 && deltaR(j, lep_m[0]) > 0.3) {
235		  if (deltaR(j, photons[0]) > 0.3 && deltaR(j, photons[1]) > 0.3 )  ++njets;
236		}
237	      }
238	      if (deltaR(j, photons[0]) < min_dR_ph0_jet) {
239		min_dR_ph0_jet = deltaR(j, photons[0]);
240		ph0_jet = j;
241	      }
242	      if (deltaR(j, photons[1]) < min_dR_ph1_jet) {
243		min_dR_ph1_jet = deltaR(j, photons[1]);
244		ph1_jet = j;
245	      }
246	    }
247	    double photon0iso = 0, photon1iso = 0;
248	    if (min_dR_ph0_jet < 0.4) photon0iso = ph0_jet.pT() - photons[0].pT();
249	    if (min_dR_ph1_jet < 0.4) photon1iso = ph1_jet.pT() - photons[1].pT();
250	    veto |= photon0iso/photons[0].pT() > 0.5;
251	    veto |= photon1iso/photons[1].pT() > 0.5;
252
253	    // Fill plots
254	    // ee and mm need doing.
255	    if (!veto) {
256	      _h["llgg"]->fill(0.5);
257	      if (el) {
258		_h["eegg"]->fill(0.5);
259	      } else {
260		_h["mmgg"]->fill(0.5);
261	      }
262
263	      if (!njets) {
264		_h["llgg"]->fill(1.5);
265		if (el) {
266		  _h["eegg"]->fill(1.5);
267		} else {
268		  _h["mmgg"]->fill(1.5);
269		}
270	      }
271	    }
272	  }
273
274	  if (deltaR(photons[0], lep_p[0]) < 0.7)  vetoEvent;
275	  if (deltaR(photons[0], lep_m[0]) < 0.7)  vetoEvent;
276
277	  // Photon isolation calculated by jets, count jets
278	  Jet ph_jet;
279	  double min_dR_ph_jet = 999.;
280	  size_t njets = 0;
281	  for (const Jet& j : jets) {
282	    if (j.pT() > 30*GeV && j.abseta() < 4.5) {
283	      if (deltaR(j, lep_p[0]) > 0.3 && deltaR(j, lep_m[0]) > 0.3 && deltaR(j, photons[0]) > 0.3)  ++njets;
284	    }
285	    if (deltaR(j, photons[0]) < min_dR_ph_jet) {
286	      min_dR_ph_jet = deltaR(j, photons[0]);
287	      ph_jet = j;
288	    }
289	  }
290
291	  double photoniso = 0;
292	  if (min_dR_ph_jet < 0.4)  photoniso = ph_jet.pT() - photons[0].pT();
293	  if (photoniso/photons[0].pT() > 0.5)  vetoEvent;
294
295
296	  // Fill plots
297	  const double pTgamma = photons[0].pT()/GeV;
298	  const double mllgamma = (lep_p[0].momentum() + lep_m[0].momentum() + photons[0].momentum()).mass()/GeV;
299
300	  _h["pT"]->fill(pTgamma);
301	  _h["M"]->fill(mllgamma);
302	  _h["Njets"]->fill(njets < 3? njets : 3);
303
304	  _h["llg"]->fill(0.5);
305	  if (el) {
306	    _h["eeg"]->fill(0.5);
307	  } else {
308	    _h["eeg"]->fill(0.5);
309	  }
310
311	  if (!njets) {
312	    _h["pT_0jet"]->fill(pTgamma);
313	    _h["M_0jet"]->fill(mllgamma);
314	    _h["llg"]->fill(1.5);
315	    if (el) {
316	      _h["eeg"]->fill(1.5);
317	    } else {
318	      _h["mmg"]->fill(1.5);
319	    }
320	  }
321	}
322      }
323    } // end of analysis
324
325
326    /// Normalise histograms etc., after the run
327    void finalize() {
328      const double sf = crossSection()/femtobarn/sumOfWeights();
329      scale(_h, sf);
330      // if we are running both e and mu, the combined lepton histos
331      // need to be divided by two to get the average
332      if (_mode == 0 || _mode == 4){
333        scale(_h["llgg"], 0.5);
334        scale(_h["llg"], 0.5);
335        scale(_h["pT"], 0.5);
336        scale(_h["pT_0jet"], 0.5);
337        scale(_h["M"], 0.5);
338        scale(_h["M_0jet"], 0.5);
339        scale(_h["Njets"], 0.5);
340      }
341    }
342
343    /// @}
344
345
346  protected:
347
348    // Data members like post-cuts event weight counters go here
349    size_t _mode;
350
351  private:
352
353    /// Histograms
354    map<string, Histo1DPtr> _h;
355
356  };
357
358
359  RIVET_DECLARE_PLUGIN(ATLAS_2016_I1448301);
360
361}