Rivet analyses referenceATLAS_2015_I1404878ttbar (to l+jets) differential cross sections at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1404878 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
Measurements of normalized differential cross-sections of top-quark pair production are presented as a function of the top-quark, $t\bar{t}$ system and event-level kinematic observables in proton-proton collisions at a centre-of-mass energy of $\sqrt{s} = 8$ TeV. The observables have been chosen to emphasize the $t\bar{t}$ production process and to be sensitive to effects of initial- and final-state radiation, to the different parton distribution functions, and to non-resonant processes and higher-order corrections. The dataset corresponds to an integrated luminosity of 20.3 fb${}^{-1}$, recorded in 2012 with the ATLAS detector at the CERN Large Hadron Collider. Events are selected in the lepton$+$jets channel, requiring exactly one charged lepton and at least four jets with at least two of the jets tagged as originating from a $b$-quark. The measured spectra are corrected for detector effects and are compared to several Monte Carlo simulations. The results are in fair agreement with the predictions over a wide kinematic range. Nevertheless, most generators predict a harder top-quark transverse momentum distribution at high values than what is observed in the data. Predictions beyond NLO accuracy improve the agreement with data at high top-quark transverse momenta. Using the current settings and parton distribution functions, the rapidity distributions are not well modelled by any generator under consideration. However, the level of agreement is improved when more recent sets of parton distribution functions are used. Source code: ATLAS_2015_I1404878.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/VetoedFinalState.hh"
5#include "Rivet/Projections/IdentifiedFinalState.hh"
6#include "Rivet/Projections/PromptFinalState.hh"
7#include "Rivet/Projections/DressedLeptons.hh"
8#include "Rivet/Projections/FastJets.hh"
9#include "Rivet/Projections/VisibleFinalState.hh"
10
11namespace Rivet {
12
13
14 class ATLAS_2015_I1404878 : public Analysis {
15 public:
16
17 /// Constructor
18 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2015_I1404878);
19
20
21 void init() {
22 // Eta ranges
23 Cut eta_full = (Cuts::abseta < 4.2) & (Cuts::pT >= 1.0*MeV);
24 Cut lep_cuts = (Cuts::abseta < 2.5) && (Cuts::pT > 25*GeV);
25
26 // All final state particles
27 FinalState fs(eta_full);
28
29 // Get photons to dress leptons
30 IdentifiedFinalState photons(fs);
31 photons.acceptIdPair(PID::PHOTON);
32
33 // Projection to find the electrons
34 IdentifiedFinalState el_id(fs);
35 el_id.acceptIdPair(PID::ELECTRON);
36
37 PromptFinalState electrons(el_id);
38 electrons.acceptTauDecays(true);
39 declare(electrons, "electrons");
40
41 DressedLeptons dressedelectrons(photons, electrons, 0.1, lep_cuts, true);
42 declare(dressedelectrons, "dressedelectrons");
43
44 DressedLeptons ewdressedelectrons(photons, electrons, 0.1, eta_full, true);
45 declare(ewdressedelectrons, "ewdressedelectrons");
46
47 // Projection to find the muons
48 IdentifiedFinalState mu_id(fs);
49 mu_id.acceptIdPair(PID::MUON);
50
51 PromptFinalState muons(mu_id);
52 muons.acceptTauDecays(true);
53 declare(muons, "muons");
54
55 DressedLeptons dressedmuons(photons, muons, 0.1, lep_cuts, true);
56 declare(dressedmuons, "dressedmuons");
57
58 DressedLeptons ewdressedmuons(photons, muons, 0.1, eta_full, true);
59 declare(ewdressedmuons, "ewdressedmuons");
60
61 // Projection to find neutrinos
62 IdentifiedFinalState nu_id;
63 nu_id.acceptNeutrinos();
64 PromptFinalState neutrinos(nu_id);
65 neutrinos.acceptTauDecays(true);
66
67 // get MET from generic invisibles
68 VetoedFinalState inv_fs(fs);
69 inv_fs.addVetoOnThisFinalState(VisibleFinalState(fs));
70 declare(inv_fs, "InvisibleFS");
71
72 // Jet clustering.
73 VetoedFinalState vfs;
74 vfs.addVetoOnThisFinalState(ewdressedelectrons);
75 vfs.addVetoOnThisFinalState(ewdressedmuons);
76 vfs.addVetoOnThisFinalState(neutrinos);
77 FastJets jets(vfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
78 declare(jets, "jets");
79
80
81 // Histogram booking
82 book(_h["massttbar"] , 1, 1, 1);
83 book(_h["massttbar_norm"] , 2, 1, 1);
84 book(_h["ptttbar"] , 3, 1, 1);
85 book(_h["ptttbar_norm"] , 4, 1, 1);
86 book(_h["absrapttbar"] , 5, 1, 1);
87 book(_h["absrapttbar_norm"] , 6, 1, 1);
88 book(_h["ptpseudotophadron"] , 7, 1, 1);
89 book(_h["ptpseudotophadron_norm"] , 8, 1, 1);
90 book(_h["absrappseudotophadron"] , 9, 1, 1);
91 book(_h["absrappseudotophadron_norm"] ,10, 1, 1);
92 book(_h["absPout"] ,11, 1, 1);
93 book(_h["absPout_norm"] ,12, 1, 1);
94 book(_h["dPhittbar"] ,13, 1, 1);
95 book(_h["dPhittbar_norm"] ,14, 1, 1);
96 book(_h["HTttbar"] ,15, 1, 1);
97 book(_h["HTttbar_norm"] ,16, 1, 1);
98 book(_h["Yboost"] ,17, 1, 1);
99 book(_h["Yboost_norm"] ,18, 1, 1);
100 book(_h["chittbar"] ,19, 1, 1);
101 book(_h["chittbar_norm"] ,20, 1, 1);
102 book(_h["RWt"] ,21, 1, 1);
103 book(_h["RWt_norm"] ,22, 1, 1);
104
105 }
106
107
108 void analyze(const Event& event) {
109
110 // Get the selected objects, using the projections.
111 vector<DressedLepton> electrons = applyProjection<DressedLeptons>(event, "dressedelectrons").dressedLeptons();
112 vector<DressedLepton> muons = applyProjection<DressedLeptons>(event, "dressedmuons").dressedLeptons();
113 const Jets& jets = applyProjection<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);
114 const FinalState& ifs = applyProjection<FinalState>(event, "InvisibleFS");
115
116 // Calculate MET
117 FourMomentum met;
118 for (const Particle& p : ifs.particles()) met += p.momentum();
119
120 // Count the number of b-tags
121 Jets bjets, lightjets;
122 for (const Jet& jet : jets){
123 bool b_tagged = jet.bTags(Cuts::pT > 5*GeV).size();
124 if ( b_tagged && bjets.size() < 2 ) bjets += jet;
125 else lightjets += jet;
126 }
127
128 bool single_electron = (electrons.size() == 1) && (muons.empty());
129 bool single_muon = (muons.size() == 1) && (electrons.empty());
130
131 DressedLepton* lepton = nullptr;
132 if (single_electron) lepton = &electrons[0];
133 else if (single_muon) lepton = &muons[0];
134
135 if(!single_electron && !single_muon) vetoEvent;
136 if (jets.size() < 4 || bjets.size() < 2) vetoEvent;
137
138 FourMomentum pbjet1; // Momentum of bjet1
139 FourMomentum pbjet2; // Momentum of bjet2
140 if ( deltaR(bjets[0], *lepton) <= deltaR(bjets[1], *lepton) ) {
141 pbjet1 = bjets[0].momentum();
142 pbjet2 = bjets[1].momentum();
143 } else {
144 pbjet1 = bjets[1].momentum();
145 pbjet2 = bjets[0].momentum();
146 }
147
148 double bestWmass = 1000.0*TeV;
149 double mWPDG = 80.399*GeV;
150 int Wj1index = -1, Wj2index = -1;
151 for (unsigned int i = 0; i < (lightjets.size() - 1); ++i) {
152 for (unsigned int j = i + 1; j < lightjets.size(); ++j) {
153 double wmass = (lightjets[i].momentum() + lightjets[j].momentum()).mass();
154 if (fabs(wmass - mWPDG) < fabs(bestWmass - mWPDG)) {
155 bestWmass = wmass;
156 Wj1index = i;
157 Wj2index = j;
158 }
159 }
160 }
161
162 // Compute hadronic W boson
163 FourMomentum pWhadron = lightjets[Wj1index].momentum() + lightjets[Wj2index].momentum();
164 double pz = _computeneutrinoz(lepton->momentum(), met);
165 FourMomentum ppseudoneutrino( sqrt(sqr(met.px()) + sqr(met.py()) + sqr(pz)), met.px(), met.py(), pz);
166
167
168 // Compute leptonic, hadronic, combined pseudo-top
169 FourMomentum ppseudotoplepton = lepton->momentum() + ppseudoneutrino + pbjet1;
170 FourMomentum ppseudotophadron = pbjet2 + pWhadron;
171 FourMomentum pttbar = ppseudotoplepton + ppseudotophadron;
172
173 Vector3 z_versor(0,0,1);
174 Vector3 vpseudotophadron = ppseudotophadron.vector3();
175 Vector3 vpseudotoplepton = ppseudotoplepton.vector3();
176 // Observables
177 double ystar = 0.5 * deltaRap(ppseudotophadron, ppseudotoplepton);
178 double chi_ttbar = exp(2 * fabs(ystar));
179 double deltaPhi_ttbar = deltaPhi(ppseudotoplepton,ppseudotophadron);
180 double HT_ttbar = ppseudotophadron.pt() + ppseudotoplepton.pt();
181 double Yboost = 0.5 * fabs(ppseudotophadron.rapidity() + ppseudotoplepton.rapidity());
182 double R_Wt = pWhadron.pt() / ppseudotophadron.pt();
183 double absPout = fabs(vpseudotophadron.dot((vpseudotoplepton.cross(z_versor))/(vpseudotoplepton.cross(z_versor).mod())));
184
185 // absolute cross sections
186 _h["ptpseudotophadron"]->fill( ppseudotophadron.pt()); //pT of pseudo top hadron
187 _h["ptttbar"]->fill( pttbar.pt()); //fill pT of ttbar in combined channel
188 _h["absrappseudotophadron"]->fill(ppseudotophadron.absrap());
189 _h["absrapttbar"]->fill( pttbar.absrap());
190 _h["massttbar"]->fill( pttbar.mass());
191 _h["absPout"]->fill( absPout);
192 _h["chittbar"]->fill( chi_ttbar);
193 _h["dPhittbar"]->fill( deltaPhi_ttbar);
194 _h["HTttbar"]->fill( HT_ttbar);
195 _h["Yboost"]->fill( Yboost);
196 _h["RWt"]->fill( R_Wt);
197 // normalised cross sections
198 _h["ptpseudotophadron_norm"]->fill( ppseudotophadron.pt()); //pT of pseudo top hadron
199 _h["ptttbar_norm"]->fill( pttbar.pt()); //fill pT of ttbar in combined channel
200 _h["absrappseudotophadron_norm"]->fill(ppseudotophadron.absrap());
201 _h["absrapttbar_norm"]->fill( pttbar.absrap());
202 _h["massttbar_norm"]->fill( pttbar.mass());
203 _h["absPout_norm"]->fill( absPout);
204 _h["chittbar_norm"]->fill( chi_ttbar);
205 _h["dPhittbar_norm"]->fill( deltaPhi_ttbar);
206 _h["HTttbar_norm"]->fill( HT_ttbar);
207 _h["Yboost_norm"]->fill( Yboost);
208 _h["RWt_norm"]->fill( R_Wt);
209
210 }
211
212
213 void finalize() {
214 // Normalize to cross-section
215 const double sf = crossSection() / sumOfWeights();
216 for (auto& k_h : _h) {
217 scale(k_h.second, sf);
218 if (k_h.first.find("_norm") != string::npos) normalize(k_h.second);
219 }
220 }
221
222
223 private:
224
225 // Compute z component of neutrino momentum given lepton and met
226 double _computeneutrinoz(const FourMomentum& lepton, FourMomentum& met) const {
227 double m_W = 80.399; // in GeV, given in the paper
228 double k = (( sqr( m_W ) - sqr( lepton.mass() ) ) / 2 ) + (lepton.px() * met.px() + lepton.py() * met.py());
229 double a = sqr ( lepton.E() )- sqr ( lepton.pz() );
230 double b = -2*k*lepton.pz();
231 double c = sqr( lepton.E() ) * sqr( met.pT() ) - sqr( k );
232 double discriminant = sqr(b) - 4 * a * c;
233 double quad[2] = { (- b - sqrt(discriminant)) / (2 * a), (- b + sqrt(discriminant)) / (2 * a) }; //two possible quadratic solns
234
235 double pzneutrino;
236 if (discriminant < 0) { // if the discriminant is negative:
237 pzneutrino = - b / (2 * a);
238 } else { // if the discriminant is positive, take the soln with smallest absolute value
239 pzneutrino = (fabs(quad[0]) < fabs(quad[1])) ? quad[0] : quad[1];
240 }
241 return pzneutrino;
242 }
243
244 /// @todo Replace with central version
245 double _mT(const FourMomentum &l, FourMomentum &nu) const {
246 return sqrt( 2 * l.pT() * nu.pT() * (1 - cos(deltaPhi(l, nu))) );
247 }
248
249
250 /// @name Objects that are used by the event selection decisions
251 map<string, Histo1DPtr> _h;
252
253 };
254
255
256 // The hook for the plugin system
257 RIVET_DECLARE_PLUGIN(ATLAS_2015_I1404878);
258
259
260}
|