Rivet analyses referenceATLAS_2012_I1203852Measurement of the $ZZ(*)$ production cross-section in $pp$ collisions at 7 TeV with ATLASExperiment: ATLAS (LHC) Inspire ID: 1203852 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Measurement of the fiducial cross section for $ZZ(*)$ production in proton proton collisions at a centre-of mass energy of 7 TeV, is presented, using data corresponding to an integrated luminosity of 4.6/fb collected by the ATLAS experiment at the Large Hadron Collider. The cross-section is measured using processes with two $Z$ bosons decaying to electrons or muons or with one $Z$ boson decaying to electrons or muons and a second $Z$ boson decaying to neutrinos. The fiducial region contains dressed leptons in restricted $p_T$ and $\eta$ ranges. The selection has specific requirements for both production processes. A measurement of the normalized fiducial cross-section as a function of $ZZ$ invariant mass, leading $Z$ $p_T$ and angle of two leptons coming from the leading $Z$ is also presented for both signal processes. Source code: ATLAS_2012_I1203852.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/IdentifiedFinalState.hh"
6#include "Rivet/Projections/PromptFinalState.hh"
7#include "Rivet/Projections/VetoedFinalState.hh"
8#include "Rivet/Projections/LeptonFinder.hh"
9#include "Rivet/Projections/MergedFinalState.hh"
10#include "Rivet/Projections/MissingMomentum.hh"
11#include "Rivet/Projections/InvMassFinalState.hh"
12
13namespace Rivet {
14
15
16 /// Generic Z candidate
17 struct Zstate : public ParticlePair {
18 Zstate() { }
19 Zstate(ParticlePair _particlepair) : ParticlePair(_particlepair) { }
20 FourMomentum mom() const { return first.momentum() + second.momentum(); }
21 operator FourMomentum() const { return mom(); }
22 static bool cmppT(const Zstate& lx, const Zstate& rx) { return lx.mom().pT() < rx.mom().pT(); }
23 };
24
25
26
27 /// ZZ analysis
28 class ATLAS_2012_I1203852 : public Analysis {
29 public:
30
31 /// Default constructor
32 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2012_I1203852);
33
34
35 void init() {
36
37 // Get options
38 // Default does everything
39 _mode = 0;
40 if ( getOption("LMODE") == "LL" ) _mode = 1;
41 if ( getOption("LMODE") == "LNU" ) _mode = 2;
42
43 // NB Missing ET is not required to be neutrinos
44 FinalState fs(Cuts::abseta < 5.0);
45 PromptFinalState pfs(fs);
46
47 // Final states to form Z bosons
48 vids.push_back(make_pair(PID::ELECTRON, PID::POSITRON));
49 vids.push_back(make_pair(PID::MUON, PID::ANTIMUON));
50
51 if (_mode!=2) {
52
53 // Selection 1: ZZ-> llll selection
54 Cut etaranges_lep = Cuts::abseta < 3.16 && Cuts::pT > 7*GeV;
55
56 LeptonFinder electron_sel4l(0.1, etaranges_lep && Cuts::abspid == PID::ELECTRON);
57 declare(electron_sel4l, "ELECTRON_sel4l");
58 LeptonFinder muon_sel4l(0.1, etaranges_lep && Cuts::abspid == PID::MUON);
59 declare(muon_sel4l, "MUON_sel4l");
60
61 // Both ZZ on-shell histos
62 book(_h_ZZ_xsect ,1, 1, 1);
63 book(_h_ZZ_ZpT ,3, 1, 1);
64 book(_h_ZZ_phill ,5, 1, 1);
65 book(_h_ZZ_mZZ ,7, 1, 1);
66
67 // One Z off-shell (ZZstar) histos
68 book(_h_ZZs_xsect ,1, 1, 2);
69
70 }
71
72 if (_mode!=1) {
73
74 // Selection 2: ZZ-> llnunu selection
75 Cut etaranges_lep2 = Cuts::abseta < 2.5 && Cuts::pT > 10*GeV;
76
77 LeptonFinder electron_sel2l2nu(0.1, etaranges_lep2 && Cuts::abspid == PID::ELECTRON);
78 declare(electron_sel2l2nu, "ELECTRON_sel2l2nu");
79 LeptonFinder muon_sel2l2nu(0.1, etaranges_lep2 && Cuts::abspid == PID::MUON);
80 declare(muon_sel2l2nu, "MUON_sel2l2nu");
81
82 /// Get all neutrinos. These will not be used to form jets.
83 IdentifiedFinalState neutrino_fs(Cuts::abseta < 4.5);
84 neutrino_fs.acceptNeutrinos();
85 declare(neutrino_fs, "NEUTRINO_FS");
86
87 // Calculate missing ET from the visible final state, not by requiring neutrinos
88 declare(MissingMomentum(Cuts::abseta < 4.5), "MISSING");
89
90 VetoedFinalState jetinput;
91 jetinput.addVetoOnThisFinalState(neutrino_fs);
92
93 FastJets jetpro(fs, JetAlg::ANTIKT, 0.4, JetMuons::NONE);
94 declare(jetpro, "jet");
95
96 // ZZ -> llnunu histos
97 book(_h_ZZnunu_xsect ,1, 1, 3);
98 book(_h_ZZnunu_ZpT ,4, 1, 1);
99 book(_h_ZZnunu_phill ,6, 1, 1);
100 book(_h_ZZnunu_mZZ ,8, 1, 1);
101
102 }
103
104 }
105
106
107 /// Do the analysis
108 void analyze(const Event& e) {
109
110 if (_mode != 2) {
111
112 ////////////////////////////////////////////////////////////////////
113 // preselection of leptons for ZZ-> llll final state
114 ////////////////////////////////////////////////////////////////////
115
116 Particles leptons_sel4l;
117
118 const DressedLeptons& mu_sel4l = apply<LeptonFinder>(e, "MUON_sel4l").dressedLeptons();
119 const DressedLeptons& el_sel4l = apply<LeptonFinder>(e, "ELECTRON_sel4l").dressedLeptons();
120
121 DressedLeptons leptonsFS_sel4l;
122 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), mu_sel4l.begin(), mu_sel4l.end() );
123 leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), el_sel4l.begin(), el_sel4l.end() );
124
125 ////////////////////////////////////////////////////////////////////
126 // OVERLAP removal dR(l,l)>0.2
127 ////////////////////////////////////////////////////////////////////
128 for (const DressedLepton& l1 : leptonsFS_sel4l) {
129 bool isolated = true;
130 for (DressedLepton& l2 : leptonsFS_sel4l) {
131 const double dR = deltaR(l1, l2);
132 if (dR < 0.2 && !isSame(l1, l2)) { isolated = false; break; }
133 }
134 if (isolated) leptons_sel4l.push_back(l1);
135 }
136
137 //////////////////////////////////////////////////////////////////
138 // Exactly two opposite charged leptons
139 //////////////////////////////////////////////////////////////////
140
141 // calculate total 'flavour' charge
142 double totalcharge = 0;
143 for (const Particle& l : leptons_sel4l) totalcharge += l.pid();
144
145 // Analyze 4 lepton events
146 if (leptons_sel4l.size() == 4 && totalcharge == 0 ) {
147 Zstate Z1, Z2;
148
149 // Identifies Z states from 4 lepton pairs
150 identifyZstates(Z1, Z2,leptons_sel4l);
151
152 ////////////////////////////////////////////////////////////////////////////
153 // Z MASS WINDOW
154 // -ZZ: for both Z: 66<mZ<116 GeV
155 // -ZZ*: one Z on-shell: 66<mZ<116 GeV, one Z off-shell: mZ>20 GeV
156 ///////////////////////////////////////////////////////////////////////////
157
158 Zstate leadPtZ = std::max(Z1, Z2, Zstate::cmppT);
159
160 double mZ1 = Z1.mom().mass();
161 double mZ2 = Z2.mom().mass();
162 double ZpT = leadPtZ.mom().pT();
163 double phill = fabs(deltaPhi(leadPtZ.first, leadPtZ.second));
164 if (phill > M_PI) phill = 2*M_PI-phill;
165 double mZZ = (Z1.mom() + Z2.mom()).mass();
166
167 if (mZ1 > 20*GeV && mZ2 > 20*GeV) {
168 // ZZ* selection
169 if (inRange(mZ1, 66*GeV, 116*GeV) || inRange(mZ2, 66*GeV, 116*GeV)) {
170 _h_ZZs_xsect -> fill(sqrtS()*GeV); ///< @todo xsec * GeV??
171 }
172
173 // ZZ selection
174 if (inRange(mZ1, 66*GeV, 116*GeV) && inRange(mZ2, 66*GeV, 116*GeV)) {
175 _h_ZZ_xsect -> fill(sqrtS()/GeV); ///< @todo xsec * GeV??
176 _h_ZZ_ZpT -> fill(ZpT);
177 _h_ZZ_phill -> fill(phill);
178 _h_ZZ_mZZ -> fill(mZZ);
179 }
180 }
181 }
182 }
183
184 if (_mode != 1) {
185
186 ////////////////////////////////////////////////////////////////////
187 /// preselection of leptons for ZZ-> llnunu final state
188 ////////////////////////////////////////////////////////////////////
189
190 Particles leptons_sel2l2nu; // output
191 const DressedLeptons& mu_sel2l2nu = apply<LeptonFinder>(e, "MUON_sel2l2nu").dressedLeptons();
192 const DressedLeptons& el_sel2l2nu = apply<LeptonFinder>(e, "ELECTRON_sel2l2nu").dressedLeptons();
193
194 DressedLeptons leptonsFS_sel2l2nu;
195 leptonsFS_sel2l2nu.insert( leptonsFS_sel2l2nu.end(), mu_sel2l2nu.begin(), mu_sel2l2nu.end() );
196 leptonsFS_sel2l2nu.insert( leptonsFS_sel2l2nu.end(), el_sel2l2nu.begin(), el_sel2l2nu.end() );
197
198 // Lepton preselection for ZZ-> llnunu
199 if ((mu_sel2l2nu.empty() || el_sel2l2nu.empty()) // cannot have opposite flavour
200 && (leptonsFS_sel2l2nu.size() == 2) // exactly two leptons
201 && (leptonsFS_sel2l2nu[0].charge() * leptonsFS_sel2l2nu[1].charge() < 1 ) // opposite charge
202 && (deltaR(leptonsFS_sel2l2nu[0], leptonsFS_sel2l2nu[1]) > 0.3) // overlap removal
203 && (leptonsFS_sel2l2nu[0].pT() > 20*GeV && leptonsFS_sel2l2nu[1].pT() > 20*GeV)) { // trigger requirement
204 leptons_sel2l2nu.insert(leptons_sel2l2nu.end(), leptonsFS_sel2l2nu.begin(), leptonsFS_sel2l2nu.end());
205 }
206 if (leptons_sel2l2nu.empty()) vetoEvent; // no further analysis, fine to veto
207
208 Particles leptons_sel2l2nu_jetveto;
209 for (const DressedLepton& l : mu_sel2l2nu) leptons_sel2l2nu_jetveto.push_back(l.bareLepton());
210 for (const DressedLepton& l : el_sel2l2nu) leptons_sel2l2nu_jetveto.push_back(l.bareLepton());
211 double ptll = (leptons_sel2l2nu[0].momentum() + leptons_sel2l2nu[1].momentum()).pT();
212
213 // Find Z1-> ll
214 FinalState fs2((Cuts::etaIn(-3.2, 3.2)));
215 InvMassFinalState imfs(fs2, vids, 20*GeV, sqrtS());
216 imfs.calc(leptons_sel2l2nu);
217 if (imfs.particlePairs().size() != 1) vetoEvent;
218 const ParticlePair& Z1constituents = imfs.particlePairs()[0];
219 FourMomentum Z1 = Z1constituents.first.momentum() + Z1constituents.second.momentum();
220
221 // Z to neutrinos candidate from missing ET
222 const MissingMomentum & missmom = apply<MissingMomentum>(e, "MISSING");
223 const FourMomentum Z2 = missmom.missingMomentum(ZMASS);
224 double met_Znunu = missmom.missingEt(); //Z2.pT();
225
226 // mTZZ
227 const double mT2_1st_term = add_quad(ZMASS, ptll) + add_quad(ZMASS, met_Znunu);
228 const double mT2_2nd_term = Z1.px() + Z2.px();
229 const double mT2_3rd_term = Z1.py() + Z2.py();
230 const double mTZZ = sqrt(sqr(mT2_1st_term) - sqr(mT2_2nd_term) - sqr(mT2_3rd_term));
231
232 if (!inRange(Z2.mass(), 66*GeV, 116*GeV)) vetoEvent;
233 if (!inRange(Z1.mass(), 76*GeV, 106*GeV)) vetoEvent;
234
235 /////////////////////////////////////////////////////////////
236 // AXIAL MET < 75 GeV
237 ////////////////////////////////////////////////////////////
238
239 double dPhiZ1Z2 = fabs(deltaPhi(Z1, Z2));
240 if (dPhiZ1Z2 > M_PI) dPhiZ1Z2 = 2*M_PI - dPhiZ1Z2;
241 const double axialEtmiss = -Z2.pT()*cos(dPhiZ1Z2);
242 if (axialEtmiss < 75*GeV) vetoEvent;
243
244 const double ZpT = Z1.pT();
245 double phill = fabs(deltaPhi(Z1constituents.first, Z1constituents.second));
246 if (phill > M_PI) phill = 2*M_PI - phill;
247
248
249 ////////////////////////////////////////////////////////////////////////////
250 // JETS
251 // -"j": found by "jetpro" projection && pT() > 25 GeV && |eta| < 4.5
252 // -"goodjets": "j" && dR(electron/muon,jet) > 0.3
253 //
254 // JETVETO: veto all events with at least one good jet
255 ///////////////////////////////////////////////////////////////////////////
256 vector<Jet> good_jets;
257 for (const Jet& j : apply<FastJets>(e, "jet").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 4.5)) {
258 bool isLepton = 0;
259 for (const Particle& l : leptons_sel2l2nu_jetveto) {
260 const double dR = deltaR(l.momentum(), j.momentum());
261 if (dR < 0.3) { isLepton = true; break; }
262 }
263 if (!isLepton) good_jets.push_back(j);
264 }
265 size_t n_sel_jets = good_jets.size();
266 if (n_sel_jets != 0) vetoEvent;
267
268
269 /////////////////////////////////////////////////////////////
270 // Fractional MET and lepton pair difference: "RatioMet"< 0.4
271 ////////////////////////////////////////////////////////////
272 double ratioMet = fabs(Z2.pT() - Z1.pT()) / Z1.pT();
273 if (ratioMet > 0.4 ) vetoEvent;
274
275
276 // End of ZZllnunu selection: now fill histograms
277 _h_ZZnunu_xsect->fill(sqrtS()/GeV); ///< @todo xsec / GeV??
278 _h_ZZnunu_ZpT ->fill(ZpT);
279 _h_ZZnunu_phill->fill(phill);
280 _h_ZZnunu_mZZ ->fill(mTZZ);
281
282 }
283 }
284
285
286 /// Finalize
287 void finalize() {
288 const double norm = crossSection()/sumOfWeights()/femtobarn;
289
290 if (_mode!=2) {
291 scale(_h_ZZ_xsect, norm);
292 normalize(_h_ZZ_ZpT);
293 normalize(_h_ZZ_phill);
294 normalize(_h_ZZ_mZZ);
295 scale(_h_ZZs_xsect, norm);
296 }
297
298 if (_mode!=1) {
299 scale(_h_ZZnunu_xsect, norm);
300 normalize(_h_ZZnunu_ZpT);
301 normalize(_h_ZZnunu_phill);
302 normalize(_h_ZZnunu_mZZ);
303 }
304 }
305
306
307 protected:
308 // Data members like post-cuts event weight counters go here
309 size_t _mode;
310
311
312 private:
313
314 void identifyZstates(Zstate& Z1, Zstate& Z2, const Particles& leptons_sel4l);
315 Histo1DPtr _h_ZZ_xsect, _h_ZZ_ZpT, _h_ZZ_phill, _h_ZZ_mZZ;
316 Histo1DPtr _h_ZZs_xsect;
317 Histo1DPtr _h_ZZnunu_xsect, _h_ZZnunu_ZpT, _h_ZZnunu_phill, _h_ZZnunu_mZZ;
318 vector< pair<PdgId,PdgId> > vids;
319 const double ZMASS = 91.1876; // GeV
320
321
322 };
323
324
325 /// 4l to ZZ assignment -- algorithm
326 void ATLAS_2012_I1203852::identifyZstates(Zstate& Z1, Zstate& Z2, const Particles& leptons_sel4l) {
327
328 /////////////////////////////////////////////////////////////////////////////
329 /// ZZ->4l pairing
330 /// - Exactly two same flavour opposite charged leptons
331 /// - Ambiguities in pairing are resolved by choosing the combination
332 /// that results in the smaller value of the sum |mll - mZ| for the two pairs
333 /////////////////////////////////////////////////////////////////////////////
334
335 Particles part_pos_el, part_neg_el, part_pos_mu, part_neg_mu;
336 for (const Particle& l : leptons_sel4l) {
337 if (l.abspid() == PID::ELECTRON) {
338 if (l.pid() < 0) part_neg_el.push_back(l);
339 if (l.pid() > 0) part_pos_el.push_back(l);
340 }
341 else if (l.abspid() == PID::MUON) {
342 if (l.pid() < 0) part_neg_mu.push_back(l);
343 if (l.pid() > 0) part_pos_mu.push_back(l);
344 }
345 }
346
347 // ee/mm channel
348 if ( part_neg_el.size() == 2 || part_neg_mu.size() == 2) {
349
350 Zstate Zcand_1, Zcand_2, Zcand_3, Zcand_4;
351 if (part_neg_el.size() == 2) { // ee
352 Zcand_1 = Zstate( ParticlePair( part_neg_el[0], part_pos_el[0] ) );
353 Zcand_2 = Zstate( ParticlePair( part_neg_el[0], part_pos_el[1] ) );
354 Zcand_3 = Zstate( ParticlePair( part_neg_el[1], part_pos_el[0] ) );
355 Zcand_4 = Zstate( ParticlePair( part_neg_el[1], part_pos_el[1] ) );
356 } else { // mumu
357 Zcand_1 = Zstate( ParticlePair( part_neg_mu[0], part_pos_mu[0] ) );
358 Zcand_2 = Zstate( ParticlePair( part_neg_mu[0], part_pos_mu[1] ) );
359 Zcand_3 = Zstate( ParticlePair( part_neg_mu[1], part_pos_mu[0] ) );
360 Zcand_4 = Zstate( ParticlePair( part_neg_mu[1], part_pos_mu[1] ) );
361 }
362
363 // We can have the following pairs: (Z1 + Z4) || (Z2 + Z3)
364 double minValue_1, minValue_2;
365 minValue_1 = fabs( Zcand_1.mom().mass() - ZMASS ) + fabs( Zcand_4.mom().mass() - ZMASS);
366 minValue_2 = fabs( Zcand_2.mom().mass() - ZMASS ) + fabs( Zcand_3.mom().mass() - ZMASS);
367 if (minValue_1 < minValue_2 ) {
368 Z1 = Zcand_1;
369 Z2 = Zcand_4;
370 } else {
371 Z1 = Zcand_2;
372 Z2 = Zcand_3;
373 }
374
375 // emu channel
376 } else if (part_neg_mu.size() == 1 && part_neg_el.size() == 1) {
377 Z1 = Zstate ( ParticlePair (part_neg_mu[0], part_pos_mu[0] ) );
378 Z2 = Zstate ( ParticlePair (part_neg_el[0], part_pos_el[0] ) );
379 }
380
381 }
382
383
384 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1203852);
385
386}
|