Rivet analyses referenceCMS_2016_I1430892Measurements of ttbar charge asymmetry using dilepton final states in pp collisions at sqrt(s) = 8 TeVExperiment: CMS (LHC) Inspire ID: 1430892 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
$\textbf{Abstract:}$ The charge asymmetry in $\mathrm{t\bar{t}}$ events is measured using dilepton final states produced in pp collisions at the LHC at $\sqrt{s}=8\:$TeV. The data sample, collected with the CMS detector, corresponds to an integrated luminosity of 19.5 $\mathrm{fb^{-1}}$. The measurements are performed using events with two oppositely charged leptons (electrons or muons) and two or more jets, where at least one of the jets is identified as originating from a bottom quark. The charge asymmetry is measured from differences in kinematic distributions, unfolded to the parton level, of positively and negatively charged top quarks and leptons. The $\mathrm{t\bar{t}}$ and leptonic charge asymmetries are found to be 0.011 +/- 0.011 (stat) +/- 0.007 (syst) and 0.003 +/- 0.006 (stat) +/- 0.003 (syst), respectively. These results, as well as charge asymmetry measurements made as a function of $\mathrm{t\bar{t}}$ system kinematic properties, are in agreement with predictions of the standard model. $\textbf{Particle-level addition to Rivet routine:}$ While the analysis was performed at the parton-level only, $\Delta|\eta_{\ell}|$ is a purely leptonic variable and it has been checked that the results of the analysis would have been essentially unchanged had it been defined at particle-level using dressed leptons instead of using the parton-level top quark daughter leptons. We therefore include both particle- and parton-level versions of this distribution in the Rivet routine, with the former identified in the plot title. For same-flavour dilepton final states, the particle-level definition in the full phase space is problematic because the two leptons can come from fully-hadronic $\mathrm{t\bar{t}}$ plus a dilepton pair from radiation. Such pairs have invariant mass $M_{\ell\ell}\sim 0$ and produce a peak near zero in the $\Delta|\eta_{\ell}|$ distribution. We therefore select only the $\mathrm{t\bar{t}}\to e\mu$ final state, by requiring exactly one electron and exactly one muon. Note this means $\mathrm{t\bar{t}}\to e\mu$ events with additional dilepton pairs from radiation are vetoed. For PYTHIA8 this amounts to 0.5% of $\mathrm{t\bar{t}}\to e\mu$ events - well below the level of sensitivity of the measured distribution. $\textbf{Histograms and covariance matrices:}$ The error bars in the measured distributions should not be used for fitting because there are significant correlations between bins. The covariance matrices for the statistical and systematic uncertainties in each distribution can be found in hepdata. The single-differential cross sections in hepdata are normalised to unit area (i.e. the integral is equal to one), while the double-differential cross sections in hepdata are normalised to the sum of entries (such that the sum of all bin heights is equal to one). This should be taken into account when comparing the measured distributions to the Rivet results and when using the covariance matrices. $\textbf{Underflow and overflow bins:}$ The lower and upper $\Delta|y_\mathrm{t}|$ and $\Delta|\eta_{\ell}|$ bins (starting and ending at -2 and +2, respectively) contain underflow and overflow events, i.e. the complete distribution from -infinity to +infinity is covered. Similarly, the upper $M_\mathrm{t\bar{t}}$, $p_\mathrm{T}^\mathrm{t\bar{t}}$, and $\left|y_\mathrm{t\bar{t}}\right|$ bins contain overflow events up to +infinity. Source code: CMS_2016_I1430892.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/IdentifiedFinalState.hh"
5#include "Rivet/Projections/PromptFinalState.hh"
6#include "Rivet/Projections/LeptonFinder.hh"
7#include "Rivet/Projections/PartonicTops.hh"
8
9namespace Rivet {
10
11
12 /// CMS 8 TeV dilepton channel ttbar charge asymmetry analysis
13 class CMS_2016_I1430892 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1430892);
18
19
20 /// Book histograms and initialise projections
21 void init() {
22
23 // Complete final state
24 FinalState fs;
25
26 // Projection for dressed electrons and muons
27 IdentifiedFinalState photons(fs);
28 photons.acceptIdPair(PID::PHOTON);
29
30 IdentifiedFinalState el_id(fs);
31 el_id.acceptIdPair(PID::ELECTRON);
32 PromptFinalState electrons(el_id);
33 declare(electrons, "Electrons");
34 LeptonFinder dressed_electrons(electrons, photons, 0.1);
35 declare(dressed_electrons, "DressedElectrons");
36
37 IdentifiedFinalState mu_id(fs);
38 mu_id.acceptIdPair(PID::MUON);
39 PromptFinalState muons(mu_id);
40 declare(muons, "Muons");
41 LeptonFinder dressed_muons(muons, photons, 0.1);
42 declare(dressed_muons, "DressedMuons");
43
44 // Parton-level top quarks
45 declare(PartonicTops(TopDecay::E_MU, PromptEMuFromTau::NO), "LeptonicPartonTops");
46
47
48 // Booking of histograms
49
50 // This histogram is independent of the parton-level information, and is an
51 // addition to the original analysis. It is compared to the same data as
52 // the parton-level delta_abseta histogram d05-x01-y01.
53 book(_h_dabsetadressedleptons, "d00-x01-y01", _bins_dabseta);
54
55 // The remaining histos use parton-level information
56 book(_h_dabseta, "d05-x01-y01", _bins_dabseta);
57 book(_h_dabsrapidity, "d02-x01-y01", _bins_dabsrapidity);
58
59 // 2D histos
60 book(_h_dabsrapidity_var[0], "d11-x01-y01", _bins_dabsrapidity, _bins_tt_mass);
61 book(_h_dabseta_var[0], "d17-x01-y01", _bins_dabseta, _bins_tt_mass);
62
63 book(_h_dabsrapidity_var[1], "d23-x01-y01", _bins_dabsrapidity, _bins_tt_pT);
64 book(_h_dabseta_var[1], "d29-x01-y01", _bins_dabseta, _bins_tt_pT);
65
66 book(_h_dabsrapidity_var[2], "d35-x01-y01", _bins_dabsrapidity, _bins_tt_absrapidity);
67 book(_h_dabseta_var[2], "d41-x01-y01", _bins_dabseta, _bins_tt_absrapidity);
68
69 // Profile histos for asymmetries
70 book(_h_dabsrapidity_profile[0], "d08-x01-y01", _bins_tt_mass);
71 book(_h_dabseta_profile[0], "d14-x01-y01", _bins_tt_mass);
72
73 book(_h_dabsrapidity_profile[1], "d20-x01-y01", _bins_tt_pT);
74 book(_h_dabseta_profile[1], "d26-x01-y01", _bins_tt_pT);
75
76 book(_h_dabsrapidity_profile[2], "d32-x01-y01", _bins_tt_absrapidity);
77 book(_h_dabseta_profile[2], "d38-x01-y01", _bins_tt_absrapidity);
78
79 }
80
81
82 /// Perform the per-event analysis
83 void analyze(const Event& event) {
84
85
86 // Use particle-level leptons for the first histogram
87 const LeptonFinder& dressed_electrons = apply<LeptonFinder>(event, "DressedElectrons");
88 const LeptonFinder& dressed_muons = apply<LeptonFinder>(event, "DressedMuons");
89
90 const DressedLeptons dressedels = dressed_electrons.dressedLeptons();
91 const DressedLeptons dressedmus = dressed_muons.dressedLeptons();
92
93 const size_t ndressedel = dressedels.size();
94 const size_t ndressedmu = dressedmus.size();
95
96 // For the particle-level histogram, require exactly one electron and exactly one muon, to select the ttbar->emu channel.
97 // Note this means ttbar->emu events with additional PromptFinalState dilepton pairs from the shower are vetoed - for PYTHIA8,
98 // this affects ~0.5% of events, so the effect is well below the level of sensitivity of the measured distribution.
99 if ( ndressedel == 1 && ndressedmu == 1 ) {
100
101 const int electrontouse = 0, muontouse = 0;
102
103 // Opposite-charge leptons only
104 if ( sameSign(dressedels[electrontouse], dressedmus[muontouse]) ) {
105 MSG_INFO("Error, e and mu have same charge, skipping event");
106 } else {
107 // Get the four-momenta of the positively- and negatively-charged leptons
108 FourMomentum lepPlus = dressedels[electrontouse].charge() > 0 ? dressedels[electrontouse] : dressedmus[muontouse];
109 FourMomentum lepMinus = dressedels[electrontouse].charge() > 0 ? dressedmus[muontouse] : dressedels[electrontouse];
110
111 // Now calculate the variable
112 double dabseta_temp = lepPlus.abseta() - lepMinus.abseta();
113
114 fillWithUFOF( _h_dabsetadressedleptons, dabseta_temp);
115 }
116
117 }
118
119
120 // The remaining variables use parton-level information.
121
122 // Get the leptonically decaying tops
123 const Particles leptonicpartontops = apply<ParticleFinder>(event, "LeptonicPartonTops").particlesByPt();
124 Particles chargedleptons;
125
126 unsigned int ntrueleptonictops = 0;
127 bool oppositesign = false;
128
129 if ( leptonicpartontops.size() == 2 ) {
130 for (size_t k = 0; k < leptonicpartontops.size(); ++k) {
131
132 // Get the lepton
133 const Particle lepTop = leptonicpartontops[k];
134 const auto isPromptChargedLepton = [](const Particle& p){return (isChargedLepton(p) && isPrompt(p, false, false));};
135 Particles lepton_candidates = lepTop.allDescendants(firstParticleWith(isPromptChargedLepton), false);
136 if ( lepton_candidates.size() < 1 ) MSG_WARNING("error, TopDecay::E_MU top quark had no daughter lepton candidate, skipping event.");
137
138 // In some cases there is no lepton from the W decay but only leptons from the decay of a radiated gamma.
139 // These hadronic PartonicTops are currently being mistakenly selected by TopDecay::E_MU (as of April 2017), and need to be rejected.
140 // TopDecay::E_MU is being fixed in Rivet, and when it is the veto below should do nothing.
141 /// @todo Should no longer be necessary -- remove
142 bool istrueleptonictop = false;
143 for (size_t i = 0; i < lepton_candidates.size(); ++i) {
144 const Particle& lepton_candidate = lepton_candidates[i];
145 if ( lepton_candidate.hasParentWith(Cuts::pid == PID::PHOTON) ) {
146 MSG_DEBUG("Found gamma parent, top: " << k+1 << " of " << leptonicpartontops.size() << " , lepton: " << i+1 << " of " << lepton_candidates.size());
147 continue;
148 }
149 if ( !istrueleptonictop && sameSign(lepTop,lepton_candidate) ) {
150 chargedleptons.push_back(lepton_candidate);
151 istrueleptonictop = true;
152 }
153 else MSG_WARNING("Error, found extra prompt charged lepton from top decay (and without gamma parent), ignoring it.");
154 }
155 if ( istrueleptonictop ) ++ntrueleptonictops;
156 }
157 }
158
159 if ( ntrueleptonictops == 2 ) {
160 oppositesign = !( sameSign(chargedleptons[0],chargedleptons[1]) );
161 if ( !oppositesign ) MSG_WARNING("error, same charge tops, skipping event.");
162 }
163
164 if ( ntrueleptonictops == 2 && oppositesign ) {
165
166 // Get the four-momenta of the positively- and negatively-charged leptons
167 const FourMomentum lepPlus = chargedleptons[0].charge() > 0 ? chargedleptons[0] : chargedleptons[1];
168 const FourMomentum lepMinus = chargedleptons[0].charge() > 0 ? chargedleptons[1] : chargedleptons[0];
169
170 const double dabseta_temp = lepPlus.abseta() - lepMinus.abseta();
171
172 // Get the four-momenta of the positively- and negatively-charged tops
173 const FourMomentum topPlus_p4 = leptonicpartontops[0].pid() > 0 ? leptonicpartontops[0] : leptonicpartontops[1];
174 const FourMomentum topMinus_p4 = leptonicpartontops[0].pid() > 0 ? leptonicpartontops[1] : leptonicpartontops[0];
175
176 const FourMomentum ttbar_p4 = topPlus_p4 + topMinus_p4;
177
178 const double tt_mass_temp = ttbar_p4.mass();
179 const double tt_absrapidity_temp = ttbar_p4.absrapidity();
180 const double tt_pT_temp = ttbar_p4.pT();
181 const double dabsrapidity_temp = topPlus_p4.absrapidity() - topMinus_p4.absrapidity();
182
183 // Fill parton-level histos
184 fillWithUFOF( _h_dabseta, dabseta_temp);
185 fillWithUFOF( _h_dabsrapidity, dabsrapidity_temp);
186
187 // Now fill the same variables in the 2D and profile histos vs ttbar invariant mass, pT, and absolute rapidity
188 for (int i_var = 0; i_var < 3; ++i_var) {
189 double var;
190 if ( i_var == 0 ) {
191 var = tt_mass_temp;
192 } else if ( i_var == 1 ) {
193 var = tt_pT_temp;
194 } else {
195 var = tt_absrapidity_temp;
196 }
197
198 fillWithUFOF( _h_dabsrapidity_var[i_var], dabsrapidity_temp, var);
199 fillWithUFOF( _h_dabseta_var[i_var], dabseta_temp, var);
200
201 fillWithUFOF( _h_dabsrapidity_profile[i_var], dabsrapidity_temp, var, (_h_dabsrapidity->xMax() + _h_dabsrapidity->xMin())/2. );
202 fillWithUFOF( _h_dabseta_profile[i_var], dabseta_temp, var, (_h_dabseta->xMax() + _h_dabseta->xMin())/2. );
203 }
204
205 }
206
207 }
208
209
210 /// Normalise histograms to unit area
211 void finalize() {
212
213 normalize(_h_dabsetadressedleptons);
214
215 normalize(_h_dabseta);
216 normalize(_h_dabsrapidity);
217
218 for (int i_var = 0; i_var < 3; ++i_var) {
219 normalize(_h_dabsrapidity_var[i_var]);
220 normalize(_h_dabseta_var[i_var]);
221 }
222
223 }
224
225
226 private:
227
228 Histo1DPtr _h_dabsetadressedleptons, _h_dabseta, _h_dabsrapidity;
229 Histo2DPtr _h_dabseta_var[3], _h_dabsrapidity_var[3];
230 Profile1DPtr _h_dabseta_profile[3], _h_dabsrapidity_profile[3];
231
232 const vector<double> _bins_tt_mass = {300., 430., 530., 1200.};
233 const vector<double> _bins_tt_pT = {0., 41., 92., 300.};
234 const vector<double> _bins_tt_absrapidity = {0., 0.34, 0.75, 1.5};
235 const vector<double> _bins_dabseta = { -2., -68./60., -48./60., -32./60., -20./60., -8./60., 0., 8./60., 20./60., 32./60., 48./60., 68./60., 2.};
236 const vector<double> _bins_dabsrapidity = {-2., -44./60., -20./60., 0., 20./60., 44./60., 2.};
237
238 void fillWithUFOF(Histo1DPtr h, double x) {
239 h->fill(std::max(std::min(x, h->xMax()-1e-9),h->xMin()+1e-9));
240 }
241
242 void fillWithUFOF(Histo2DPtr h, double x, double y) {
243 h->fill(std::max(std::min(x, h->xMax()-1e-9),h->xMin()+1e-9), std::max(std::min(y, h->yMax()-1e-9),h->yMin()+1e-9));
244 }
245
246 void fillWithUFOF(Profile1DPtr h, double x, double y, double c) {
247 h->fill(std::max(std::min(y, h->xMax()-1e-9),h->xMin()+1e-9), float(x > c) - float(x < c));
248 }
249
250
251 };
252
253
254 RIVET_DECLARE_PLUGIN(CMS_2016_I1430892);
255
256
257}
|