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/LeptonFinder.hh"
8#include "Rivet/Projections/FastJets.hh"
9#include "Rivet/Projections/VisibleFinalState.hh"
10
11namespace Rivet {
12
13
14 /// ttbar (to l+jets) differential cross sections at 8 TeV
15 class ATLAS_2015_I1404878 : public Analysis {
16 public:
17
18 /// Constructor
19 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2015_I1404878);
20
21
22 void init() {
23 // Eta ranges
24 Cut eta_full = (Cuts::abseta < 4.2) & (Cuts::pT >= 1.0*MeV);
25 Cut lep_cuts = (Cuts::abseta < 2.5) && (Cuts::pT > 25*GeV);
26
27 // All final state particles
28 FinalState fs(eta_full);
29
30 // Get photons to dress leptons
31 IdentifiedFinalState photons(fs);
32 photons.acceptIdPair(PID::PHOTON);
33
34 // Projection to find the electrons
35 IdentifiedFinalState el_id(fs);
36 el_id.acceptIdPair(PID::ELECTRON);
37
38 PromptFinalState electrons(el_id);
39 electrons.acceptTauDecays(true);
40 declare(electrons, "electrons");
41
42 LeptonFinder dressedelectrons(electrons, photons, 0.1, lep_cuts);
43 declare(dressedelectrons, "dressedelectrons");
44
45 LeptonFinder ewdressedelectrons(electrons, photons, 0.1, eta_full);
46 declare(ewdressedelectrons, "ewdressedelectrons");
47
48 // Projection to find the muons
49 IdentifiedFinalState mu_id(fs);
50 mu_id.acceptIdPair(PID::MUON);
51
52 PromptFinalState muons(mu_id);
53 muons.acceptTauDecays(true);
54 declare(muons, "muons");
55
56 LeptonFinder dressedmuons(muons, photons, 0.1, lep_cuts);
57 declare(dressedmuons, "dressedmuons");
58
59 LeptonFinder ewdressedmuons(muons, photons, 0.1, eta_full);
60 declare(ewdressedmuons, "ewdressedmuons");
61
62 // Projection to find neutrinos
63 IdentifiedFinalState nu_id;
64 nu_id.acceptNeutrinos();
65 PromptFinalState neutrinos(nu_id);
66 neutrinos.acceptTauDecays(true);
67
68 // get MET from generic invisibles
69 VetoedFinalState inv_fs(fs);
70 inv_fs.addVetoOnThisFinalState(VisibleFinalState(fs));
71 declare(inv_fs, "InvisibleFS");
72
73 // Jet clustering.
74 VetoedFinalState vfs;
75 vfs.addVetoOnThisFinalState(ewdressedelectrons);
76 vfs.addVetoOnThisFinalState(ewdressedmuons);
77 vfs.addVetoOnThisFinalState(neutrinos);
78 FastJets jets(vfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::DECAY);
79 declare(jets, "jets");
80
81
82 // Histogram booking
83 book(_h["massttbar"] , 1, 1, 1);
84 book(_h["massttbar_norm"] , 2, 1, 1);
85 book(_h["ptttbar"] , 3, 1, 1);
86 book(_h["ptttbar_norm"] , 4, 1, 1);
87 book(_h["absrapttbar"] , 5, 1, 1);
88 book(_h["absrapttbar_norm"] , 6, 1, 1);
89 book(_h["ptpseudotophadron"] , 7, 1, 1);
90 book(_h["ptpseudotophadron_norm"] , 8, 1, 1);
91 book(_h["absrappseudotophadron"] , 9, 1, 1);
92 book(_h["absrappseudotophadron_norm"] ,10, 1, 1);
93 book(_h["absPout"] ,11, 1, 1);
94 book(_h["absPout_norm"] ,12, 1, 1);
95 book(_h["dPhittbar"] ,13, 1, 1);
96 book(_h["dPhittbar_norm"] ,14, 1, 1);
97 book(_h["HTttbar"] ,15, 1, 1);
98 book(_h["HTttbar_norm"] ,16, 1, 1);
99 book(_h["Yboost"] ,17, 1, 1);
100 book(_h["Yboost_norm"] ,18, 1, 1);
101 book(_h["chittbar"] ,19, 1, 1);
102 book(_h["chittbar_norm"] ,20, 1, 1);
103 book(_h["RWt"] ,21, 1, 1);
104 book(_h["RWt_norm"] ,22, 1, 1);
105
106 }
107
108
109 void analyze(const Event& event) {
110
111 // Get the selected objects, using the projections.
112 DressedLeptons electrons = apply<LeptonFinder>(event, "dressedelectrons").dressedLeptons();
113 DressedLeptons muons = apply<LeptonFinder>(event, "dressedmuons").dressedLeptons();
114 const Jets& jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);
115 const FinalState& ifs = apply<FinalState>(event, "InvisibleFS");
116
117 // Calculate MET
118 FourMomentum met;
119 for (const Particle& p : ifs.particles()) met += p.momentum();
120
121 // Count the number of b-tags
122 Jets bjets, lightjets;
123 for (const Jet& jet : jets){
124 bool b_tagged = jet.bTags(Cuts::pT > 5*GeV).size();
125 if ( b_tagged && bjets.size() < 2 ) bjets += jet;
126 else lightjets += jet;
127 }
128
129 bool single_electron = (electrons.size() == 1) && (muons.empty());
130 bool single_muon = (muons.size() == 1) && (electrons.empty());
131
132 DressedLepton* lepton = nullptr;
133 if (single_electron) lepton = &electrons[0];
134 else if (single_muon) lepton = &muons[0];
135
136 if(!single_electron && !single_muon) vetoEvent;
137 if (jets.size() < 4 || bjets.size() < 2) vetoEvent;
138
139 FourMomentum pbjet1; // Momentum of bjet1
140 FourMomentum pbjet2; // Momentum of bjet2
141 if ( deltaR(bjets[0], *lepton) <= deltaR(bjets[1], *lepton) ) {
142 pbjet1 = bjets[0].momentum();
143 pbjet2 = bjets[1].momentum();
144 } else {
145 pbjet1 = bjets[1].momentum();
146 pbjet2 = bjets[0].momentum();
147 }
148
149 double bestWmass = 1000.0*TeV;
150 double mWPDG = 80.399*GeV;
151 int Wj1index = -1, Wj2index = -1;
152 for (unsigned int i = 0; i < (lightjets.size() - 1); ++i) {
153 for (unsigned int j = i + 1; j < lightjets.size(); ++j) {
154 double wmass = (lightjets[i].momentum() + lightjets[j].momentum()).mass();
155 if (fabs(wmass - mWPDG) < fabs(bestWmass - mWPDG)) {
156 bestWmass = wmass;
157 Wj1index = i;
158 Wj2index = j;
159 }
160 }
161 }
162
163 // Compute hadronic W boson
164 FourMomentum pWhadron = lightjets[Wj1index].momentum() + lightjets[Wj2index].momentum();
165 double pz = _computeneutrinoz(lepton->momentum(), met);
166 FourMomentum ppseudoneutrino( sqrt(sqr(met.px()) + sqr(met.py()) + sqr(pz)), met.px(), met.py(), pz);
167
168
169 // Compute leptonic, hadronic, combined pseudo-top
170 FourMomentum ppseudotoplepton = lepton->momentum() + ppseudoneutrino + pbjet1;
171 FourMomentum ppseudotophadron = pbjet2 + pWhadron;
172 FourMomentum pttbar = ppseudotoplepton + ppseudotophadron;
173
174 Vector3 z_versor(0,0,1);
175 Vector3 vpseudotophadron = ppseudotophadron.vector3();
176 Vector3 vpseudotoplepton = ppseudotoplepton.vector3();
177 // Observables
178 double ystar = 0.5 * deltaRap(ppseudotophadron, ppseudotoplepton);
179 double chi_ttbar = exp(2 * fabs(ystar));
180 double deltaPhi_ttbar = deltaPhi(ppseudotoplepton,ppseudotophadron);
181 double HT_ttbar = ppseudotophadron.pt() + ppseudotoplepton.pt();
182 double Yboost = 0.5 * fabs(ppseudotophadron.rapidity() + ppseudotoplepton.rapidity());
183 double R_Wt = pWhadron.pt() / ppseudotophadron.pt();
184 double absPout = fabs(vpseudotophadron.dot((vpseudotoplepton.cross(z_versor))/(vpseudotoplepton.cross(z_versor).mod())));
185
186 // absolute cross sections
187 _h["ptpseudotophadron"]->fill( ppseudotophadron.pt()); //pT of pseudo top hadron
188 _h["ptttbar"]->fill( pttbar.pt()); //fill pT of ttbar in combined channel
189 _h["absrappseudotophadron"]->fill(ppseudotophadron.absrap());
190 _h["absrapttbar"]->fill( pttbar.absrap());
191 _h["massttbar"]->fill( pttbar.mass());
192 _h["absPout"]->fill( absPout);
193 _h["chittbar"]->fill( chi_ttbar);
194 _h["dPhittbar"]->fill( deltaPhi_ttbar);
195 _h["HTttbar"]->fill( HT_ttbar);
196 _h["Yboost"]->fill( Yboost);
197 _h["RWt"]->fill( R_Wt);
198 // normalised cross sections
199 _h["ptpseudotophadron_norm"]->fill( ppseudotophadron.pt()); //pT of pseudo top hadron
200 _h["ptttbar_norm"]->fill( pttbar.pt()); //fill pT of ttbar in combined channel
201 _h["absrappseudotophadron_norm"]->fill(ppseudotophadron.absrap());
202 _h["absrapttbar_norm"]->fill( pttbar.absrap());
203 _h["massttbar_norm"]->fill( pttbar.mass());
204 _h["absPout_norm"]->fill( absPout);
205 _h["chittbar_norm"]->fill( chi_ttbar);
206 _h["dPhittbar_norm"]->fill( deltaPhi_ttbar);
207 _h["HTttbar_norm"]->fill( HT_ttbar);
208 _h["Yboost_norm"]->fill( Yboost);
209 _h["RWt_norm"]->fill( R_Wt);
210
211 }
212
213
214 void finalize() {
215 // Normalize to cross-section
216 const double sf = crossSection()/picobarn / sumOfWeights();
217 for (auto& k_h : _h) {
218 scale(k_h.second, sf);
219 if (k_h.first.find("_norm") != string::npos) normalize(k_h.second);
220 }
221 }
222
223
224 private:
225
226 // Compute z component of neutrino momentum given lepton and met
227 double _computeneutrinoz(const FourMomentum& lepton, FourMomentum& met) const {
228 double m_W = 80.399; // in GeV, given in the paper
229 double k = (( sqr( m_W ) - sqr( lepton.mass() ) ) / 2 ) + (lepton.px() * met.px() + lepton.py() * met.py());
230 double a = sqr ( lepton.E() )- sqr ( lepton.pz() );
231 double b = -2*k*lepton.pz();
232 double c = sqr( lepton.E() ) * sqr( met.pT() ) - sqr( k );
233 double discriminant = sqr(b) - 4 * a * c;
234 double quad[2] = { (- b - sqrt(discriminant)) / (2 * a), (- b + sqrt(discriminant)) / (2 * a) }; //two possible quadratic solns
235
236 double pzneutrino;
237 if (discriminant < 0) { // if the discriminant is negative:
238 pzneutrino = - b / (2 * a);
239 } else { // if the discriminant is positive, take the soln with smallest absolute value
240 pzneutrino = (fabs(quad[0]) < fabs(quad[1])) ? quad[0] : quad[1];
241 }
242 return pzneutrino;
243 }
244
245 /// @todo Replace with central version
246 double _mT(const FourMomentum &l, FourMomentum &nu) const {
247 return sqrt( 2 * l.pT() * nu.pT() * (1 - cos(deltaPhi(l, nu))) );
248 }
249
250
251 /// @name Objects that are used by the event selection decisions
252 map<string, Histo1DPtr> _h;
253
254 };
255
256
257 RIVET_DECLARE_PLUGIN(ATLAS_2015_I1404878);
258
259
260}
|