Rivet analyses referenceATLAS_2014_I1306615Higgs diphoton events at 8 TeV in ATLASExperiment: ATLAS (LHC) Inspire ID: 1306615 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
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}
|