rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1306615

Higgs diphoton events at 8 TeV in ATLAS
Experiment: ATLAS (LHC)
Inspire ID: 1306615
Status: VALIDATED
Authors:
  • Michaela Queitsch-Maitland
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • pp Higgs production at 8 TeV, include all processes (ggH, VH, VBF, ttH)

Measurements of fiducial and differential cross sections are presented for Higgs boson production in proton-proton collisions at a centre-of-mass energy of 8TeV. The analysis is performed in the Higgs diphoton decay channel using 20.3/fb of data recorded by the ATLAS experiment at the CERN Large Hadron Collider. The signal yields are corrected for the effects of detector inefficiency and resolution. Differential cross sections are also presented, as a function of variables related to the diphoton kinematics and the jet activity produced in the Higgs boson events. The observed spectra are statistically limited but broadly in line with the theoretical expectations. Note on reinterpretation: Since the background is subtracted from a fit, any BSM signal would presumably also have beed subtracted. Therefore should not be used to constrain BSM signals unless they comes from a genuine excess of SM Higgs->gammagamma.

Source code: ATLAS_2014_I1306615.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/PromptFinalState.hh"
  5#include "Rivet/Projections/LeptonFinder.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7#include "Rivet/Projections/FastJets.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief ATLAS H->yy differential cross-sections measurement
 13  ///
 14  /// @author Michaela Queitsch-Maitland <michaela.queitsch-maitland@cern.ch>
 15  //
 16  // arXiv: http://arxiv.org/abs/ARXIV:1407.4222
 17  // HepData: http://hepdata.cedar.ac.uk/view/ins1306615
 18  class ATLAS_2014_I1306615 : public Analysis {
 19  public:
 20
 21    // Constructor
 22    ATLAS_2014_I1306615()
 23      : Analysis("ATLAS_2014_I1306615")
 24    {    }
 25
 26
 27    // Book histograms and initialise projections before the run
 28    void init() {
 29
 30      // Final state
 31      // All particles within |eta| < 5.0
 32      const FinalState FS(Cuts::abseta<5.0);
 33      declare(FS,"FS");
 34
 35      // Project photons with pT > 25 GeV and |eta| < 2.37
 36      PromptFinalState ph_FS(Cuts::abseta<2.37 && Cuts::pT>25*GeV && Cuts::pid == PID::PHOTON);
 37      declare(ph_FS, "PH_FS");
 38
 39      // Project photons for dressing
 40      FinalState ph_dressing_FS(Cuts::abspid == PID::PHOTON);
 41
 42      // Project bare electrons
 43      PromptFinalState el_bare_FS(Cuts::abseta < 5.0 && Cuts::abspid == PID::ELECTRON);
 44
 45      // Project dressed electrons with pT > 15 GeV and |eta| < 2.47
 46      LeptonFinder el_dressed_FS(el_bare_FS, ph_dressing_FS, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 15*GeV);
 47      declare(el_dressed_FS,"EL_DRESSED_FS");
 48
 49      // Project bare muons
 50      PromptFinalState mu_bare_FS(Cuts::abseta < 5.0 && Cuts::abspid == PID::MUON);
 51
 52      // Project dressed muons with pT > 15 GeV and |eta| < 2.47
 53      LeptonFinder mu_dressed_FS(mu_bare_FS, ph_dressing_FS, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 15*GeV);
 54      declare(mu_dressed_FS,"MU_DRESSED_FS");
 55
 56      // Final state excluding muons and neutrinos (for jet building and photon isolation)
 57      VetoedFinalState veto_mu_nu_FS(FS);
 58      veto_mu_nu_FS.vetoNeutrinos();
 59      veto_mu_nu_FS.addVetoPairId(PID::MUON);
 60      declare(veto_mu_nu_FS, "VETO_MU_NU_FS");
 61
 62      // Build the anti-kT R=0.4 jets, using FinalState particles (vetoing muons and neutrinos)
 63      FastJets jets(veto_mu_nu_FS, JetAlg::ANTIKT, 0.4);
 64      declare(jets, "JETS");
 65
 66      // Book histograms
 67      // 1D distributions
 68      book(_h_pT_yy         ,1,1,1);
 69      book(_h_y_yy          ,2,1,1);
 70      book(_h_Njets30       ,3,1,1);
 71      book(_h_Njets50       ,4,1,1);
 72      book(_h_pT_j1         ,5,1,1);
 73      book(_h_y_j1          ,6,1,1);
 74      book(_h_HT            ,7,1,1);
 75      book(_h_pT_j2         ,8,1,1);
 76      book(_h_Dy_jj         ,9,1,1);
 77      book(_h_Dphi_yy_jj    ,10,1,1);
 78      book(_h_cosTS_CS      ,11,1,1);
 79      book(_h_cosTS_CS_5bin ,12,1,1);
 80      book(_h_Dphi_jj       ,13,1,1);
 81      book(_h_pTt_yy        ,14,1,1);
 82      book(_h_Dy_yy         ,15,1,1);
 83      book(_h_tau_jet       ,16,1,1);
 84      book(_h_sum_tau_jet   ,17,1,1);
 85      book(_h_y_j2          ,18,1,1);
 86      book(_h_pT_j3         ,19,1,1);
 87      book(_h_m_jj          ,20,1,1);
 88      book(_h_pT_yy_jj      ,21,1,1);
 89
 90      // 2D distributions of cosTS_CS x pT_yy
 91      book(_h_cosTS_pTyy_low,  22,1,1);
 92      book(_h_cosTS_pTyy_high, 22,1,2);
 93      book(_h_cosTS_pTyy_rest, 22,1,3);
 94
 95      // 2D distributions of Njets x pT_yy
 96      book(_h_pTyy_Njets0, 23,1,1);
 97      book(_h_pTyy_Njets1, 23,1,2);
 98      book(_h_pTyy_Njets2, 23,1,3);
 99
100      book(_h_pTj1_excl, 24,1,1);
101
102      // Fiducial regions
103      book(_h_fidXSecs, 30,1,1);
104    }
105
106    // Perform the per-event analysis
107    void analyze(const Event& event) {
108
109      // Get final state particles
110      const Particles& FS_ptcls         = apply<FinalState>(event, "FS").particles();
111      const Particles& ptcls_veto_mu_nu = apply<VetoedFinalState>(event, "VETO_MU_NU_FS").particles();
112      const Particles& photons          = apply<PromptFinalState>(event, "PH_FS").particlesByPt();
113      DressedLeptons good_el          = apply<LeptonFinder>(event, "EL_DRESSED_FS").dressedLeptons();
114      DressedLeptons good_mu          = apply<LeptonFinder>(event, "MU_DRESSED_FS").dressedLeptons();
115
116      // For isolation calculation
117      float dR_iso    = 0.4;
118      float ETcut_iso = 14.0;
119      FourMomentum ET_iso;
120
121      // Fiducial selection: pT > 25 GeV, |eta| < 2.37 and isolation (in cone deltaR = 0.4) is < 14 GeV
122      Particles fid_photons;
123      for (const Particle& ph : photons) {
124
125        // Calculate isolation
126        ET_iso = - ph.momentum();
127        // Loop over fs truth particles (excluding muons and neutrinos)
128        for (const Particle& p : ptcls_veto_mu_nu) {
129          // Check if the truth particle is in a cone of 0.4
130          if ( deltaR(ph, p) < dR_iso ) ET_iso += p.momentum();
131        }
132
133        // Check isolation
134        if ( ET_iso.Et() > ETcut_iso ) continue;
135
136        // Fill vector of photons passing fiducial selection
137        fid_photons.push_back(ph);
138      }
139
140      if (fid_photons.size() < 2)  vetoEvent;
141
142      const FourMomentum y1 = fid_photons[0].momentum();
143      const FourMomentum y2 = fid_photons[1].momentum();
144
145      double m_yy = (y1 + y2).mass();
146
147      // Relative pT cuts
148      if ( y1.pT() < 0.35 * m_yy || y2.pT() < 0.25 * m_yy ) vetoEvent;
149
150      // Mass window cut
151      if ( m_yy < 105 || m_yy > 160 ) vetoEvent;
152
153      // -------------------------------------------- //
154      // Passed diphoton baseline fiducial selection! //
155      // -------------------------------------------- //
156
157      // Muon and Electron selection
158      idiscard(good_mu, [&](const DressedLepton& lep) { return deltaR(lep, y1) < 0.4 || deltaR(lep, y2) < 0.4; });
159      idiscard(good_el, [&](const DressedLepton& lep) { return deltaR(lep, y1) < 0.4 || deltaR(lep, y2) < 0.4; });
160
161      // Find prompt, invisible particles for missing ET calculation
162      // Based on VisibleFinalState projection
163      FourMomentum invisible(0,0,0,0);
164      for (const Particle& p : FS_ptcls) {
165
166        // Veto non-prompt particles (from hadron or tau decay)
167        if ( !p.isPrompt() ) continue;
168        // Charged particles are visible
169        if ( PID::charge3( p.pid() ) != 0 ) continue;
170        // Neutral hadrons are visible
171        if ( PID::isHadron( p.pid() ) ) continue;
172        // Photons are visible
173        if ( p.pid() == PID::PHOTON ) continue;
174        // Gluons are visible (for parton level analyses)
175        if ( p.pid() == PID::GLUON ) continue;
176        // Everything else is invisible
177        invisible += p.momentum();
178      }
179      double MET = invisible.Et();
180
181      // Jet selection
182      // Get jets with pT > 25 GeV and |rapidity| < 4.4
183      const Jets& jets = apply<FastJets>(event, "JETS").jetsByPt(Cuts::pT>25*GeV && Cuts::absrap <4.4);
184
185      Jets jets_25, jets_30, jets_50;
186
187      for (const Jet& jet : jets) {
188
189        bool passOverlap = true;
190        // Overlap with leading photons
191        if ( deltaR(y1, jet.momentum()) < 0.4 ) passOverlap = false;
192        if ( deltaR(y2, jet.momentum()) < 0.4 ) passOverlap = false;
193
194        // Overlap with good electrons
195        for (const auto& el : good_el) {
196          if ( deltaR(el, jet) < 0.2 ) passOverlap = false;
197        }
198
199        if ( ! passOverlap ) continue;
200
201        if ( jet.abseta() < 2.4 || ( jet.abseta() > 2.4 && jet.pT() > 30*GeV) ) jets_25 += jet;
202        if ( jet.pT() > 30*GeV ) jets_30 += jet;
203        if ( jet.pT() > 50*GeV ) jets_50 += jet;
204      }
205
206      // Fiducial regions
207      _h_fidXSecs->fill(1);
208      if ( jets_30.size() >= 1 ) _h_fidXSecs->fill(2);
209      if ( jets_30.size() >= 2 ) _h_fidXSecs->fill(3);
210      if ( jets_30.size() >= 3 ) _h_fidXSecs->fill(4);
211      if ( jets_30.size() >= 2 && passVBFCuts(y1 + y2, jets_30.at(0).momentum(), jets_30.at(1).momentum()) ) _h_fidXSecs->fill(5);
212      if ( (good_el.size() + good_mu.size()) > 0 ) _h_fidXSecs->fill(6);
213      if ( MET > 80 ) _h_fidXSecs->fill(7);
214
215      // Fill histograms
216      // Inclusive variables
217      _pT_yy    = (y1 + y2).pT();
218      _y_yy     = (y1 + y2).absrap();
219      _cosTS_CS = cosTS_CS(y1, y2);
220      _pTt_yy   = pTt(y1, y2);
221      _Dy_yy    = fabs( deltaRap(y1, y2) );
222
223      _Njets30 = jets_30.size() > 3 ? 3 : jets_30.size();
224      _Njets50 = jets_50.size() > 3 ? 3 : jets_50.size();
225      _h_Njets30->fill(_Njets30);
226      _h_Njets50->fill(_Njets50);
227
228      _pT_j1 = jets_30.size() > 0 ? jets_30[0].momentum().pT() : 0.;
229      _pT_j2 = jets_30.size() > 1 ? jets_30[1].momentum().pT() : 0.;
230      _pT_j3 = jets_30.size() > 2 ? jets_30[2].momentum().pT() : 0.;
231
232      _HT = 0.0;
233      for (const Jet& jet : jets_30) {  _HT += jet.pT(); }
234
235      _tau_jet     = tau_jet_max(y1 + y2, jets_25);
236      _sum_tau_jet = sum_tau_jet(y1 + y2, jets_25);
237
238      _h_pT_yy        ->fill(_pT_yy);
239      _h_y_yy         ->fill(_y_yy);
240      _h_pT_j1        ->fill(_pT_j1);
241      _h_cosTS_CS     ->fill(_cosTS_CS);
242      _h_cosTS_CS_5bin->fill(_cosTS_CS);
243      _h_HT           ->fill(_HT);
244      _h_pTt_yy       ->fill(_pTt_yy);
245      _h_Dy_yy        ->fill(_Dy_yy);
246      _h_tau_jet      ->fill(_tau_jet);
247      _h_sum_tau_jet  ->fill(_sum_tau_jet);
248
249      // >=1 jet variables
250      if ( jets_30.size() >= 1 ) {
251        FourMomentum j1 = jets_30[0].momentum();
252        _y_j1 = j1.absrap();
253
254	_h_pT_j2->fill(_pT_j2);
255	_h_y_j1 ->fill(_y_j1);
256      }
257
258      // >=2 jet variables
259      if ( jets_30.size() >= 2 ) {
260        FourMomentum j1 = jets_30[0].momentum();
261        FourMomentum j2 = jets_30[1].momentum();
262
263        _Dy_jj      = fabs( deltaRap(j1, j2) );
264        _Dphi_jj    = fabs( deltaPhi(j1, j2) );
265        _Dphi_yy_jj = fabs( deltaPhi(y1 + y2, j1 + j2) );
266        _m_jj       = (j1 + j2).mass();
267        _pT_yy_jj   = (y1 + y2 + j1 + j2).pT();
268        _y_j2       = j2.absrap();
269
270	_h_Dy_jj      ->fill(_Dy_jj);
271	_h_Dphi_jj    ->fill(_Dphi_jj);
272	_h_Dphi_yy_jj ->fill(_Dphi_yy_jj);
273	_h_m_jj       ->fill(_m_jj);
274	_h_pT_yy_jj   ->fill(_pT_yy_jj);
275	_h_pT_j3      ->fill(_pT_j3);
276	_h_y_j2       ->fill(_y_j2);
277      }
278
279      // 2D distributions of cosTS_CS x pT_yy
280      if ( _pT_yy < 80 )
281	_h_cosTS_pTyy_low->fill(_cosTS_CS);
282      else if ( _pT_yy > 80 && _pT_yy < 200 )
283	_h_cosTS_pTyy_high->fill(_cosTS_CS);
284      else if ( _pT_yy > 200 )
285	_h_cosTS_pTyy_rest->fill(_cosTS_CS);
286
287      // 2D distributions of pT_yy x Njets
288      if ( _Njets30 == 0 )
289	_h_pTyy_Njets0->fill(_pT_yy);
290      else if ( _Njets30 == 1 )
291	_h_pTyy_Njets1->fill(_pT_yy);
292      else if ( _Njets30 >= 2 )
293	_h_pTyy_Njets2->fill(_pT_yy);
294
295      if ( _Njets30 == 1 ) _h_pTj1_excl->fill(_pT_j1);
296
297    }
298
299    // Normalise histograms after the run
300    void finalize() {
301
302      const double xs = crossSectionPerEvent()/femtobarn;
303
304      scale(_h_pT_yy, xs);
305      scale(_h_y_yy, xs);
306      scale(_h_pT_j1, xs);
307      scale(_h_y_j1, xs);
308      scale(_h_HT, xs);
309      scale(_h_pT_j2, xs);
310      scale(_h_Dy_jj, xs);
311      scale(_h_Dphi_yy_jj, xs);
312      scale(_h_cosTS_CS, xs);
313      scale(_h_cosTS_CS_5bin, xs);
314      scale(_h_Dphi_jj, xs);
315      scale(_h_pTt_yy, xs);
316      scale(_h_Dy_yy, xs);
317      scale(_h_tau_jet, xs);
318      scale(_h_sum_tau_jet, xs);
319      scale(_h_y_j2, xs);
320      scale(_h_pT_j3, xs);
321      scale(_h_m_jj, xs);
322      scale(_h_pT_yy_jj, xs);
323      scale(_h_cosTS_pTyy_low, xs);
324      scale(_h_cosTS_pTyy_high, xs);
325      scale(_h_cosTS_pTyy_rest, xs);
326      scale(_h_pTyy_Njets0, xs);
327      scale(_h_pTyy_Njets1, xs);
328      scale(_h_pTyy_Njets2, xs);
329      scale(_h_pTj1_excl, xs);
330      scale(_h_Njets30, xs);
331      scale(_h_Njets50, xs);
332      scale(_h_fidXSecs, xs);
333    }
334
335    // VBF-enhanced dijet topology selection cuts
336    bool passVBFCuts(const FourMomentum &H, const FourMomentum &j1, const FourMomentum &j2) {
337      return ( fabs(deltaRap(j1, j2)) > 2.8 && (j1 + j2).mass() > 400 && fabs(deltaPhi(H, j1 + j2)) > 2.6 );
338    }
339
340    // Cosine of the decay angle in the Collins-Soper frame
341    double cosTS_CS(const FourMomentum &y1, const FourMomentum &y2) {
342      return fabs( ( (y1.E() + y1.pz())* (y2.E() - y2.pz()) - (y1.E() - y1.pz()) * (y2.E() + y2.pz()) )
343		   / ((y1 + y2).mass() * sqrt(pow((y1 + y2).mass(), 2) + pow((y1 + y2).pt(), 2)) ) );
344    }
345
346    // Diphoton pT along thrust axis
347    double pTt(const FourMomentum &y1, const FourMomentum &y2) {
348      return fabs(y1.px() * y2.py() - y2.px() * y1.py()) / (y1 - y2).pT()*2;
349    }
350
351    // Tau of jet  (see paper for description)
352    // tau_jet = mT/(2*cosh(y*)), where mT = pT (+) m, and y* = rapidty in Higgs rest frame
353    double tau_jet( const FourMomentum &H, const FourMomentum &jet ) {
354      return sqrt( pow(jet.pT(),2) + pow(jet.mass(),2) ) / (2.0 * cosh( jet.rapidity() - H.rapidity() ) );
355    }
356
357    // Maximal (leading) tau_jet (see paper for description)
358    double tau_jet_max(const FourMomentum &H, const Jets& jets, double tau_jet_cut = 8.) {
359      double max_tj = 0;
360      for (const auto& jet : jets) {
361        FourMomentum j = jet.momentum();
362        if (tau_jet(H, j) > tau_jet_cut)  max_tj = max(tau_jet(H, j), max_tj);
363      }
364      return max_tj;
365    }
366
367    // Scalar sum of tau for all jets (see paper for description)
368    double sum_tau_jet(const FourMomentum &H, const Jets& jets, double tau_jet_cut = 8.)  {
369      double sum_tj = 0;
370      for (const auto& jet : jets) {
371        FourMomentum j = jet.momentum();
372        if (tau_jet(H, j) > tau_jet_cut)  sum_tj += tau_jet(H, j);
373      }
374      return sum_tj;
375    }
376
377  private:
378
379    Histo1DPtr _h_pT_yy;
380    Histo1DPtr _h_y_yy;
381    Histo1DPtr _h_Njets30;
382    Histo1DPtr _h_Njets50;
383    Histo1DPtr _h_pT_j1;
384    Histo1DPtr _h_y_j1;
385    Histo1DPtr _h_HT;
386    Histo1DPtr _h_pT_j2;
387    Histo1DPtr _h_Dy_jj;
388    Histo1DPtr _h_Dphi_yy_jj;
389    Histo1DPtr _h_cosTS_CS;
390    Histo1DPtr _h_cosTS_CS_5bin;
391    Histo1DPtr _h_Dphi_jj;
392    Histo1DPtr _h_pTt_yy;
393    Histo1DPtr _h_Dy_yy;
394    Histo1DPtr _h_tau_jet;
395    Histo1DPtr _h_sum_tau_jet;
396    Histo1DPtr _h_y_j2;
397    Histo1DPtr _h_pT_j3;
398    Histo1DPtr _h_m_jj;
399    Histo1DPtr _h_pT_yy_jj;
400    Histo1DPtr _h_cosTS_pTyy_low;
401    Histo1DPtr _h_cosTS_pTyy_high;
402    Histo1DPtr _h_cosTS_pTyy_rest;
403    Histo1DPtr _h_pTyy_Njets0;
404    Histo1DPtr _h_pTyy_Njets1;
405    Histo1DPtr _h_pTyy_Njets2;
406    Histo1DPtr _h_pTj1_excl;
407    Histo1DPtr _h_fidXSecs;
408
409    int _Njets30;
410    int _Njets50;
411    double _pT_yy;
412    double _y_yy;
413    double _cosTS_CS;
414    double _pT_j1;
415    double _m_jj;
416    double _y_j1;
417    double _HT;
418    double _pT_j2;
419    double _y_j2;
420    double _Dphi_yy_jj;
421    double _pT_yy_jj;
422    double _Dphi_jj;
423    double _Dy_jj;
424    double _pT_j3;
425    double _pTt_yy;
426    double _Dy_yy;
427    double _tau_jet;
428    double _sum_tau_jet;
429  };
430
431
432  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1306615);
433
434}