Rivet analyses referenceATLAS_2024_I2771257Z+b(b) and Z+c at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 2771257 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
This paper presents a measurement of the production cross-section of a $Z$ boson in association with $b$- or $c$-jets, in proton–proton collisions at $\sqrt{s} =$13 TeV with the ATLAS experiment at the Large Hadron Collider using data corresponding to an integrated luminosity of 140 fb$^{-1}$. Inclusive and differential cross-sections are measured for events containing a $Z$ boson decaying into electrons or muons and produced in association with at least one $b$-jet, at least one $c$-jet, or at least two $b$-jets with transverse momentum $p_\text{T} > 20$GeV and rapidity $|y| < 2.5$. Predictions from several Monte Carlo generators based on next-to-leading-order matrix elements interfaced with a parton-shower simulation, with different choices of flavour schemes for initial-state partons, are compared with the measured cross-sections. The results are also compared with novel predictions, based on infrared and collinear safe jet flavour dressing algorithms. Selected $Z+\geq 1$ $c$-jet observables, optimized for sensitivity to intrinsic-charm, are compared with benchmark models with different intrinsic-charm fractions. Source code: ATLAS_2024_I2771257.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/DileptonFinder.hh"
4#include "Rivet/Projections/VetoedFinalState.hh"
5#include "Rivet/Projections/PromptFinalState.hh"
6#include "Rivet/Projections/HeavyHadrons.hh"
7#include "Rivet/Projections/FastJets.hh"
8
9namespace Rivet {
10
11 /// @brief Z+b(b) and Z+c at 13 TeV
12 class ATLAS_2024_I2771257 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2024_I2771257);
17
18 /// Book histograms and initialise projections before the run
19 void init() {
20
21 // Initialise and register projections
22
23 _mode = 0;
24 if ( getOption("LMODE") == "EL" ) _mode = 1;
25 if ( getOption("LMODE") == "MU" ) _mode = 2;
26
27 const FinalState fs;
28 // Define fiducial cuts for the leptons in the ZFinder
29 const Cut lepcuts = (Cuts::pT > 27*GeV) && (Cuts::abseta < 2.5);
30 DileptonFinder zfinderE(91.2*GeV, 0.1, lepcuts && Cuts::abspid == PID::ELECTRON, Cuts::massIn(76*GeV, 106*GeV));
31 DileptonFinder zfinderM(91.2*GeV, 0.1, lepcuts && Cuts::abspid == PID::MUON, Cuts::massIn(76*GeV, 106*GeV));
32 declare(zfinderE, "zfinderE");
33 declare(zfinderM, "zfinderM");
34
35 //Build jets using ATLAS AntiKt4WZtruthJets definition removing dressed W, Z leptons
36
37 // Photons
38 FinalState photons(Cuts::abspid == PID::PHOTON);
39
40 // Muons
41 PromptFinalState bare_mu(Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT);
42 LeptonFinder all_dressed_mu(bare_mu, photons, 0.1, Cuts::abseta < 2.5);
43
44 // Electrons
45 PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT);
46 LeptonFinder all_dressed_el(bare_el, photons, 0.1, Cuts::abseta < 2.5);
47
48 //Jet forming
49 VetoedFinalState vfs(FinalState(Cuts::abseta < 4.5));
50 vfs.addVetoOnThisFinalState(all_dressed_el);
51 vfs.addVetoOnThisFinalState(all_dressed_mu);
52
53 FastJets jetfs(vfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::DECAY);
54 declare(jetfs, "jets");
55
56 //HF hadrons for b and c jet def
57 declare(HeavyHadrons(), "HFHadrons");
58
59 // fiducial1B
60 book(_h["fiducial1B_Z_Pt"], 4, 1, 1);
61 book(_h["fiducial1B_leadBJet_Pt"], 5, 1, 1);
62 book(_h["fiducial1B_ZleadBJet_DR"], 6, 1, 1);
63
64 // fiducial2B
65 book(_h["fiducial2B_diBJets_DPhi"], 7, 1, 1);
66 book(_h["fiducial2B_diBJets_M"], 8, 1, 1);
67
68 // fiducial1C
69 book(_h["fiducial1C_Z_Pt"], 9, 1, 1);
70 book(_h["fiducial1C_leadCJet_Pt"], 10, 1, 1);
71 book(_h["fiducial1C_leadCJet_xF"], 11, 1, 1);
72
73 }
74
75 /// Perform the per-event analysis
76 void analyze(const Event& event) {
77
78 // --------------------------------------- object sel: leptons -----------------------------------------------
79
80 const DileptonFinder& zfinderE = apply<DileptonFinder>(event, "zfinderE");
81 const Particles& els = zfinderE.constituents();
82 const DileptonFinder& zfinderM = apply<DileptonFinder>(event, "zfinderM");
83 const Particles& mus = zfinderM.constituents();
84
85 // default is to run average of Z->ee and Z->mm
86 // use LMODE option to pick one channel
87 if ( (els.size() + mus.size()) != 2 ) vetoEvent;
88
89 if ( _mode == 0 && !(els.size()==2 || mus.size()==2) ) vetoEvent;
90 else if ( _mode == 1 && !(els.size() == 2 && mus.empty()) ) vetoEvent;
91 else if ( _mode == 2 && !(els.empty() && mus.size() == 2) ) vetoEvent;
92
93 double Z_pT = 0., Zrap = 0., Zphi = 0.;
94 if ( els.size()==2 ) {
95 Z_pT = zfinderE.boson().pt()/GeV;
96 Zphi = zfinderE.boson().phi();
97 Zrap = zfinderE.boson().rapidity();
98 } else {
99 Z_pT = zfinderM.boson().pt()/GeV;
100 Zphi = zfinderM.boson().phi();
101 Zrap = zfinderM.boson().rapidity();
102 }
103
104 // --------------------------------------- object sel: jets ---------------------------------------------------
105 // Retrieve clustered jets, sorted by pT, with a minimum pT cut
106 Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 20*GeV && Cuts::absrap < 2.5);
107
108 // --------------------------------------- object sel: OR ---------------------------------------------------
109 //OR 0.4 between jets and dressed leptons
110 idiscardIfAnyDeltaRLess(jets, els, 0.4);
111 idiscardIfAnyDeltaRLess(jets, mus, 0.4);
112
113 //b- and c- jet definition according to cone R=0.4 labelling
114 Jets bjets, cjets;
115 const HeavyHadrons &ha = apply<HeavyHadrons>(event,"HFHadrons");
116 const Particles allBs = ha.bHadrons(Cuts::pT > 5*GeV);
117 const Particles allCs = ha.cHadrons(Cuts::pT > 5*GeV);
118 Particles matchedBs, matchedCs;
119
120 std::vector<int> bJet_idx_v;
121 int jet_idx =0;
122
123 for (const Jet& j : jets) {
124 Jet closest_j;
125 Particle closest_b;
126 double minDR_j_b = 10;
127
128 for (const Particle& bHad : allBs) {
129 bool alreadyMatched = false;
130 for (const Particle& bMatched : matchedBs) {
131 alreadyMatched |= bMatched.isSame(bHad);
132 }
133 if (alreadyMatched) continue;
134
135 double DR_j_b = deltaR(j, bHad);
136 if (DR_j_b <= 0.3 && DR_j_b < minDR_j_b) {
137 minDR_j_b = DR_j_b;
138 closest_j = j;
139 closest_b = bHad;
140 }
141 }
142
143 if (minDR_j_b < 0.3) {
144 bjets += closest_j;
145 matchedBs += closest_b;
146 bJet_idx_v.push_back(jet_idx);
147 }
148 ++jet_idx;
149 }
150
151 jet_idx = 0;
152 for (const Jet& j : jets) {
153 Jet closest_j;
154 Particle closest_c;
155 double minDR_j_c = 10;
156
157 for (const Particle& cHad : allCs) {
158 bool alreadyMatched = false;
159 for (const Particle& cMatched : matchedCs) {
160 alreadyMatched |= cMatched.isSame(cHad);
161 }
162
163 // Do not match c-jets if already b-matched
164 for (size_t bJet_i = 0; bJet_i<bJet_idx_v.size(); ++bJet_i) {
165 if (jet_idx == bJet_idx_v[bJet_i]) alreadyMatched = true;
166 }
167
168 if (alreadyMatched) continue;
169
170 double DR_j_c = deltaR(j, cHad);
171 if ( DR_j_c <= 0.3 && DR_j_c < minDR_j_c) {
172 minDR_j_c = DR_j_c;
173 closest_j = j;
174 closest_c = cHad;
175 }
176 }
177
178 if (minDR_j_c < 0.3) {
179 cjets += closest_j;
180 matchedCs += closest_c;
181 }
182 ++jet_idx;
183 }
184
185 // --------------------------------------- event sel: Nbjets ------------------------------------------------------
186 if (bjets.size() < 1 && cjets.size() < 1) vetoEvent;
187
188 // --------------------------------------- event classification ------------------------------------------------------
189 // string region = "";
190 bool pass_fid1B = false;
191 bool pass_fid1C = false;
192 bool pass_fid2B = false;
193
194 // based on the flavour of the leading HF jet
195 if ( (bjets.size()>=1 && cjets.size()==0) ||
196 (bjets.size()>=1 && cjets.size()>=1 && bjets[0].pT() > cjets[0].pT()) ) {
197 pass_fid1B = true;
198 }
199 if ( (cjets.size()>=1 && bjets.size()==0) ||
200 (cjets.size()>=1 && bjets.size()>=1 && bjets[0].pT() < cjets[0].pT()) ) {
201 pass_fid1C = true;
202 }
203
204 // based on the flavour of the leading HF jet
205 if ( (bjets.size()>=2 && cjets.size()==0) ||
206 (bjets.size()>=2 && cjets.size()>=1 && bjets[0].pT() > cjets[0].pT()) ) {
207 pass_fid2B = true;
208 }
209 // -------------------------------------------------------------------------------------------------------------------
210
211
212
213 // --------------------------------------- observable cal -------------------------------------------------------
214
215 double leadBJet_Pt = -100, ZleadBJet_DR = -100, dRapZb = -100, dPhiZb = -100, diBJets_M = -100, diBJets_DPhi = -100;
216
217 if (pass_fid1B){
218 leadBJet_Pt = bjets[0].pT();
219 dRapZb=fabs(Zrap-bjets[0].rap());
220 dPhiZb=acos(cos(fabs(Zphi-bjets[0].phi())));
221 ZleadBJet_DR = sqrt(dRapZb*dRapZb+dPhiZb*dPhiZb);
222 }
223
224 if (pass_fid2B){
225 diBJets_M = (bjets[0].momentum() + bjets[1].momentum()).mass();
226 diBJets_DPhi = acos(cos(fabs(bjets[0].phi()-bjets[1].phi())));
227 }
228
229 // Feynam scaling variable defined as (2 pT*sinh(eta))/sqrt(s) will have the negative range
230 // and pT*sinh(eta) == pZ, which is not exactly accurate for massive objects, such as jet
231 double leadCJet_Pt = -100, leadCJet_xF = -100;
232
233 if (pass_fid1C){
234 leadCJet_Pt = cjets[0].pT();
235 leadCJet_xF = 2*fabs(cjets[0].pz()/GeV)/13000.;
236 }
237
238 // ------------------------------------------ hist fill ---------------------------------------------------------
239 if (pass_fid1B){
240 _h["fiducial1B_Z_Pt"]->fill(Z_pT/GeV);
241 _h["fiducial1B_leadBJet_Pt"]->fill(leadBJet_Pt/GeV);
242 _h["fiducial1B_ZleadBJet_DR"]->fill(ZleadBJet_DR);
243 }
244
245 if (pass_fid2B){
246 _h["fiducial2B_diBJets_M"]->fill(diBJets_M/GeV);
247 _h["fiducial2B_diBJets_DPhi"]->fill(diBJets_DPhi);
248 }
249
250 if (pass_fid1C){
251 _h["fiducial1C_Z_Pt"]->fill(Z_pT/GeV);
252 _h["fiducial1C_leadCJet_Pt"]->fill(leadCJet_Pt/GeV);
253 _h["fiducial1C_leadCJet_xF"]->fill(leadCJet_xF);
254 }
255
256 // -------------------------------------------------------------------------------------------------------------
257
258 }
259
260
261 /// Normalise histograms etc., after the run
262 void finalize() {
263 const double sf = _mode? 1.0 : 0.5;
264 scale(_h, sf * crossSectionPerEvent());
265 }
266
267 private:
268
269 size_t _mode;
270
271 map<string, Histo1DPtr> _h;
272
273 };
274
275
276 RIVET_DECLARE_PLUGIN(ATLAS_2024_I2771257);
277
278}
|