Rivet analyses referenceATLAS_2014_I1310835H(125) -> 4l at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1310835 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
Measurements of fiducial and differential cross sections of Higgs boson production in the $H\to ZZ^\ast\to 4\ell$ decay channel are presented. The cross sections are determined within a fiducial phase space and corrected for detection efficiency and resolution effects. They are based on 20.3 fb$^{-1}$ of $pp$ collision data, produced at $\sqrt{s}=8$ TeV centre-of-mass energy at the LHC and recorded by the ATLAS detector. The differential measurements are performed in bins of transverse momentum and rapidity of the four-lepton system, the invariant mass of the subleading lepton pair and the decay angle of the leading lepton pair with respect to the beam line in the four-lepton rest frame, as well as the number of jets and the transverse momentum of the leading jet. The measured cross sections are compared to selected theoretical calculations of the Standard Model expectations. No significant deviation from any of the tested predictions is found. Source code: ATLAS_2014_I1310835.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/PromptFinalState.hh"
6#include "Rivet/Projections/LeptonFinder.hh"
7
8namespace Rivet {
9
10
11 /// @brief H(125)->ZZ->4l at 8 TeV
12 class ATLAS_2014_I1310835 : public Analysis {
13 public:
14
15 /// Default constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1310835);
17
18 void init() {
19 const FinalState fs(Cuts::abseta < 5.0);
20
21 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
22
23 PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON);
24
25 PromptFinalState bare_mu(Cuts::abspid == PID::MUON);
26
27 // Selection: lepton selection
28 Cut etaranges_el = Cuts::abseta < 2.47 && Cuts::pT > 7*GeV;
29 LeptonFinder electron_sel4l(bare_el, photons, 0.1, etaranges_el);
30 declare(electron_sel4l, "electrons");
31
32 Cut etaranges_mu = Cuts::abseta < 2.7 && Cuts::pT > 6*GeV;
33 LeptonFinder muon_sel4l(bare_mu, photons, 0.1, etaranges_mu);
34 declare(muon_sel4l, "muons");
35
36 FastJets jetpro(fs, JetAlg::ANTIKT, 0.4, JetMuons::NONE, JetInvisibles::NONE);
37 declare(jetpro, "jet");
38
39 // Book histos
40 book(_h_pt , 1, 1, 1);
41 book(_h_rapidity , 2, 1, 1);
42 book(_h_m34 , 3, 1, 1);
43 book(_h_costheta , 4, 1, 1);
44 book(_h_njets , 5, 1, 1);
45 book(_h_leadingjetpt, 6, 1, 1);
46
47 }
48
49
50
51 /// Do the analysis
52 void analyze(const Event& e) {
53
54 if (_sedges.empty()) _sedges = _h_njets->xEdges();
55 ////////////////////////////////////////////////////////////////////
56 // preselection of leptons for ZZ-> llll final state
57 ////////////////////////////////////////////////////////////////////
58
59 const DressedLeptons& mu_sel4l = apply<LeptonFinder>(e, "muons").dressedLeptons();
60 const DressedLeptons& el_sel4l = apply<LeptonFinder>(e, "electrons").dressedLeptons();
61
62 DressedLeptons leptonsFS_sel4l;
63 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), mu_sel4l.begin(), mu_sel4l.end() );
64 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), el_sel4l.begin(), el_sel4l.end() );
65
66 /////////////////////////////////////////////////////////////////////////////
67 /// H->ZZ->4l pairing
68 /////////////////////////////////////////////////////////////////////////////
69
70 size_t el_p = 0;
71 size_t el_n = 0;
72 size_t mu_p = 0;
73 size_t mu_n = 0;
74
75 for (const Particle& l : leptonsFS_sel4l) {
76 if (l.abspid() == PID::ELECTRON) {
77 if (l.pid() < 0) ++el_n;
78 if (l.pid() > 0) ++el_p;
79 }
80 else if (l.abspid() == PID::MUON) {
81 if (l.pid() < 0) ++mu_n;
82 if (l.pid() > 0) ++mu_p;
83 }
84 }
85
86 bool pass_sfos = ( (el_p >=2 && el_n >=2) || (mu_p >=2 && mu_n >=2) || (el_p >=1 && el_n >=1 && mu_p >=1 && mu_n >=1) );
87
88 if (!pass_sfos) vetoEvent;
89
90 Zstate Z1, Z2, Zcand;
91 size_t n_parts = leptonsFS_sel4l.size();
92 size_t l1_index = 0;
93 size_t l2_index = 0;
94
95 // determine Z1 first
96 double min_mass_diff = -1;
97 for (size_t i = 0; i < n_parts; ++i) {
98 for (size_t j = 0; j < n_parts; ++j) {
99 if (i >= j) continue;
100
101 if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid()) continue; //only pair SFOS leptons
102
103 Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
104 double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
105
106 if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
107 min_mass_diff = mass_diff;
108 Z1 = Zcand;
109 l1_index = i;
110 l2_index = j;
111 }
112 }
113 }
114
115 //determine Z2 second
116 min_mass_diff = -1;
117 for (size_t i = 0; i < n_parts; ++i) {
118 if (i == l1_index || i == l2_index) continue;
119 for (size_t j = 0; j < n_parts; ++j) {
120 if (j == l1_index || j == l2_index || i >= j) continue;
121
122 if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid()) continue; // only pair SFOS leptons
123
124 Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
125 double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
126
127 if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
128 min_mass_diff = mass_diff;
129 Z2 = Zcand;
130 }
131 }
132 }
133
134 Particles leptons_sel4l;
135 leptons_sel4l.push_back(Z1.first);
136 leptons_sel4l.push_back(Z1.second);
137 leptons_sel4l.push_back(Z2.first);
138 leptons_sel4l.push_back(Z2.second);
139
140 ////////////////////////////////////////////////////////////////////////////
141 // Kinematic Requirements
142 ///////////////////////////////////////////////////////////////////////////
143
144 //leading lepton pT requirement
145 std::vector<double> lepton_pt;
146 for (const Particle& i : leptons_sel4l) lepton_pt.push_back(i.pT() / GeV);
147 std::sort(lepton_pt.begin(), lepton_pt.end(), [](const double pT1, const double pT2) { return pT1 > pT2; });
148
149 if (!(lepton_pt[0] > 20*GeV && lepton_pt[1] > 15*GeV && lepton_pt[2] > 10*GeV)) vetoEvent;
150
151 //invariant mass requirements
152 if (!(inRange(Z1.mom().mass(), 50*GeV, 106*GeV) && inRange(Z2.mom().mass(), 12*GeV, 115*GeV))) vetoEvent;
153
154 //lepton separation requirements
155 for (unsigned int i = 0; i < 4; ++i) {
156 for (unsigned int j = 0; j < 4; ++j) {
157 if (i >= j) continue;
158 double dR = deltaR(leptons_sel4l[i], leptons_sel4l[j]);
159 bool sameflavor = leptons_sel4l[i].abspid() == leptons_sel4l[j].abspid();
160
161 if ( sameflavor && dR < 0.1) vetoEvent;
162 if (!sameflavor && dR < 0.2) vetoEvent;
163 }
164 }
165
166 // J/Psi veto requirement
167 for (unsigned int i = 0; i < 4; ++i) {
168 for (unsigned int j = 0; j < 4; ++j) {
169 if (i >= j) continue;
170 if ( leptons_sel4l[i].pid() != -1*leptons_sel4l[j].pid() ) continue;
171 if ((leptons_sel4l[i].momentum() + leptons_sel4l[j].momentum()).mass() <= 5*GeV) vetoEvent;
172 }
173 }
174
175 // 4-lepton invariant mass requirement
176 double m4l = (Z1.mom() + Z2.mom()).mass();
177 if (!(inRange(m4l, 118*GeV, 129*GeV))) vetoEvent;
178
179
180 ////////////////////////////////////////////////////////////////////////////
181 // Higgs observables
182 ///////////////////////////////////////////////////////////////////////////
183 FourMomentum Higgs = Z1.mom() + Z2.mom();
184
185 double H4l_pt = Higgs.pt()/GeV;
186 double H4l_rapidity = Higgs.absrap();
187 LorentzTransform HRF_boost;
188 //HRF_boost.mkFrameTransformFromBeta(Higgs.betaVec());
189 HRF_boost.setBetaVec(- Higgs.betaVec());
190 FourMomentum Z1_in_HRF = HRF_boost.transform( Z1.mom() );
191 double H4l_costheta = fabs(cos( Z1_in_HRF.theta()));
192 double H4l_m34 = Z2.mom().mass()/GeV;
193
194 ////////////////////////////////////////////////////////////////////////////
195 // Jet observables
196 ///////////////////////////////////////////////////////////////////////////
197
198 Jets jets;
199 for (const Jet& jet : apply<FastJets>(e, "jet").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.4)) {
200 bool overlaps = false;
201 for (const Particle& lep : leptonsFS_sel4l) {
202 if (lep.abspid() != PID::ELECTRON) continue;
203 const double dR = deltaR(lep, jet);
204 if (dR < 0.2) { overlaps = true; break; }
205 }
206 if (!overlaps) jets += jet;
207 }
208 size_t n_jets = jets.size();
209 if (n_jets > 3) n_jets = 3;
210
211 std::vector<double> jet_pt;
212 for (const Jet& i : jets) jet_pt.push_back(i.pT()/GeV);
213
214 double leading_jet_pt = n_jets? jet_pt[0] : 0.;
215
216 ////////////////////////////////////////////////////////////////////////////
217 // End of H->ZZ->llll selection: now fill histograms
218 ////////////////////////////////////////////////////////////////////////////
219
220
221 _h_pt->fill(H4l_pt);
222 _h_rapidity->fill(H4l_rapidity);
223 _h_costheta->fill(H4l_costheta);
224 _h_m34->fill(H4l_m34);
225 _h_njets->fill(_sedges[min(3, n_jets)]);
226 _h_leadingjetpt->fill(leading_jet_pt);
227
228
229 }
230
231
232 /// Generic Z candidate
233 struct Zstate : public ParticlePair {
234 Zstate() { }
235 Zstate(ParticlePair _particlepair) : ParticlePair(_particlepair) { }
236 FourMomentum mom() const { return first.momentum() + second.momentum(); }
237 operator FourMomentum() const { return mom(); }
238 };
239
240 /// Finalize
241 void finalize() {
242
243 const double norm = crossSection()/sumOfWeights()/femtobarn;
244 scale(_h_pt, norm);
245 scale(_h_rapidity, norm);
246 scale(_h_costheta, norm);
247 scale(_h_m34, norm);
248 scale(_h_njets, norm);
249 scale(_h_leadingjetpt, norm);
250 }
251
252
253 private:
254
255 Histo1DPtr _h_pt, _h_rapidity, _h_costheta;
256 Histo1DPtr _h_m34, _h_leadingjetpt;
257 BinnedHistoPtr<string> _h_njets;
258 vector<string> _sedges;
259
260 };
261
262
263 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1310835);
264
265}
|