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 ////////////////////////////////////////////////////////////////////
55 // preselection of leptons for ZZ-> llll final state
56 ////////////////////////////////////////////////////////////////////
57
58 const DressedLeptons& mu_sel4l = apply<LeptonFinder>(e, "muons").dressedLeptons();
59 const DressedLeptons& el_sel4l = apply<LeptonFinder>(e, "electrons").dressedLeptons();
60
61 DressedLeptons leptonsFS_sel4l;
62 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), mu_sel4l.begin(), mu_sel4l.end() );
63 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), el_sel4l.begin(), el_sel4l.end() );
64
65 /////////////////////////////////////////////////////////////////////////////
66 /// H->ZZ->4l pairing
67 /////////////////////////////////////////////////////////////////////////////
68
69 size_t el_p = 0;
70 size_t el_n = 0;
71 size_t mu_p = 0;
72 size_t mu_n = 0;
73
74 for (const Particle& l : leptonsFS_sel4l) {
75 if (l.abspid() == PID::ELECTRON) {
76 if (l.pid() < 0) ++el_n;
77 if (l.pid() > 0) ++el_p;
78 }
79 else if (l.abspid() == PID::MUON) {
80 if (l.pid() < 0) ++mu_n;
81 if (l.pid() > 0) ++mu_p;
82 }
83 }
84
85 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) );
86
87 if (!pass_sfos) vetoEvent;
88
89 Zstate Z1, Z2, Zcand;
90 size_t n_parts = leptonsFS_sel4l.size();
91 size_t l1_index = 0;
92 size_t l2_index = 0;
93
94 // determine Z1 first
95 double min_mass_diff = -1;
96 for (size_t i = 0; i < n_parts; ++i) {
97 for (size_t j = 0; j < n_parts; ++j) {
98 if (i >= j) continue;
99
100 if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid()) continue; //only pair SFOS leptons
101
102 Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
103 double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
104
105 if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
106 min_mass_diff = mass_diff;
107 Z1 = Zcand;
108 l1_index = i;
109 l2_index = j;
110 }
111 }
112 }
113
114 //determine Z2 second
115 min_mass_diff = -1;
116 for (size_t i = 0; i < n_parts; ++i) {
117 if (i == l1_index || i == l2_index) continue;
118 for (size_t j = 0; j < n_parts; ++j) {
119 if (j == l1_index || j == l2_index || i >= j) continue;
120
121 if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid()) continue; // only pair SFOS leptons
122
123 Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
124 double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
125
126 if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
127 min_mass_diff = mass_diff;
128 Z2 = Zcand;
129 }
130 }
131 }
132
133 Particles leptons_sel4l;
134 leptons_sel4l.push_back(Z1.first);
135 leptons_sel4l.push_back(Z1.second);
136 leptons_sel4l.push_back(Z2.first);
137 leptons_sel4l.push_back(Z2.second);
138
139 ////////////////////////////////////////////////////////////////////////////
140 // Kinematic Requirements
141 ///////////////////////////////////////////////////////////////////////////
142
143 //leading lepton pT requirement
144 std::vector<double> lepton_pt;
145 for (const Particle& i : leptons_sel4l) lepton_pt.push_back(i.pT() / GeV);
146 std::sort(lepton_pt.begin(), lepton_pt.end(), [](const double pT1, const double pT2) { return pT1 > pT2; });
147
148 if (!(lepton_pt[0] > 20*GeV && lepton_pt[1] > 15*GeV && lepton_pt[2] > 10*GeV)) vetoEvent;
149
150 //invariant mass requirements
151 if (!(inRange(Z1.mom().mass(), 50*GeV, 106*GeV) && inRange(Z2.mom().mass(), 12*GeV, 115*GeV))) vetoEvent;
152
153 //lepton separation requirements
154 for (unsigned int i = 0; i < 4; ++i) {
155 for (unsigned int j = 0; j < 4; ++j) {
156 if (i >= j) continue;
157 double dR = deltaR(leptons_sel4l[i], leptons_sel4l[j]);
158 bool sameflavor = leptons_sel4l[i].abspid() == leptons_sel4l[j].abspid();
159
160 if ( sameflavor && dR < 0.1) vetoEvent;
161 if (!sameflavor && dR < 0.2) vetoEvent;
162 }
163 }
164
165 // J/Psi veto requirement
166 for (unsigned int i = 0; i < 4; ++i) {
167 for (unsigned int j = 0; j < 4; ++j) {
168 if (i >= j) continue;
169 if ( leptons_sel4l[i].pid() != -1*leptons_sel4l[j].pid() ) continue;
170 if ((leptons_sel4l[i].momentum() + leptons_sel4l[j].momentum()).mass() <= 5*GeV) vetoEvent;
171 }
172 }
173
174 // 4-lepton invariant mass requirement
175 double m4l = (Z1.mom() + Z2.mom()).mass();
176 if (!(inRange(m4l, 118*GeV, 129*GeV))) vetoEvent;
177
178
179 ////////////////////////////////////////////////////////////////////////////
180 // Higgs observables
181 ///////////////////////////////////////////////////////////////////////////
182 FourMomentum Higgs = Z1.mom() + Z2.mom();
183
184 double H4l_pt = Higgs.pt()/GeV;
185 double H4l_rapidity = Higgs.absrap();
186 LorentzTransform HRF_boost;
187 //HRF_boost.mkFrameTransformFromBeta(Higgs.betaVec());
188 HRF_boost.setBetaVec(- Higgs.betaVec());
189 FourMomentum Z1_in_HRF = HRF_boost.transform( Z1.mom() );
190 double H4l_costheta = fabs(cos( Z1_in_HRF.theta()));
191 double H4l_m34 = Z2.mom().mass()/GeV;
192
193 ////////////////////////////////////////////////////////////////////////////
194 // Jet observables
195 ///////////////////////////////////////////////////////////////////////////
196
197 Jets jets;
198 for (const Jet& jet : apply<FastJets>(e, "jet").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.4)) {
199 bool overlaps = false;
200 for (const Particle& lep : leptonsFS_sel4l) {
201 if (lep.abspid() != PID::ELECTRON) continue;
202 const double dR = deltaR(lep, jet);
203 if (dR < 0.2) { overlaps = true; break; }
204 }
205 if (!overlaps) jets += jet;
206 }
207 size_t n_jets = jets.size();
208 if (n_jets > 3) n_jets = 3;
209
210 std::vector<double> jet_pt;
211 for (const Jet& i : jets) jet_pt.push_back(i.pT()/GeV);
212
213 double leading_jet_pt = n_jets? jet_pt[0] : 0.;
214
215 ////////////////////////////////////////////////////////////////////////////
216 // End of H->ZZ->llll selection: now fill histograms
217 ////////////////////////////////////////////////////////////////////////////
218
219
220 _h_pt->fill(H4l_pt);
221 _h_rapidity->fill(H4l_rapidity);
222 _h_costheta->fill(H4l_costheta);
223 _h_m34->fill(H4l_m34);
224 _h_njets->fill(n_jets + 1);
225 _h_leadingjetpt->fill(leading_jet_pt);
226
227
228 }
229
230
231 /// Generic Z candidate
232 struct Zstate : public ParticlePair {
233 Zstate() { }
234 Zstate(ParticlePair _particlepair) : ParticlePair(_particlepair) { }
235 FourMomentum mom() const { return first.momentum() + second.momentum(); }
236 operator FourMomentum() const { return mom(); }
237 };
238
239 /// Finalize
240 void finalize() {
241
242 const double norm = crossSection()/sumOfWeights()/femtobarn;
243 scale(_h_pt, norm);
244 scale(_h_rapidity, norm);
245 scale(_h_costheta, norm);
246 scale(_h_m34, norm);
247 scale(_h_njets, norm);
248 scale(_h_leadingjetpt, norm);
249 }
250
251
252 private:
253
254 Histo1DPtr _h_pt, _h_rapidity, _h_costheta;
255 Histo1DPtr _h_m34, _h_njets, _h_leadingjetpt;
256
257 };
258
259
260 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1310835);
261
262}
|