rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2015_I1345452

Pseudo-top-antitop cross sections
Experiment: ATLAS (LHC)
Inspire ID: 1345452
Status: VALIDATED
Authors:
  • Longen Lan
  • Aldo Saavedra
  • Kevin Finelli
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • top-antitop production, reconstructed from leptons + jets final state

Various differential cross-sections are measured in top-quark pair ($t\bar{t}$) events produced in proton-proton collisions at a centre-of-mass energy of $\sqrt{s}=7$ TeV at the LHC with the ATLAS detector. These differential cross-sections are presented in a data set corresponding to an integrated luminosity of 4.6 fb$^{-1}$. The differential cross-sections are presented in terms of kinematic variables of a top-quark proxy referred to as the pseudo-top-quark whose dependence on theoretical models is minimal. The pseudo-top-quark can be defined in terms of either reconstructed detector objects or stable particles in an analogous way. The measurements are performed on $t\bar{t}$ events in the lepton+jets channel, requiring exactly one charged lepton and at least four jets with at least two of them tagged as originating from a $b$-quark. The hadronic and leptonic pseudo-top-quarks are defined via the leptonic or hadronic decay mode of the W boson produced by the top-quark decay in events with a single charged lepton. The cross-section is measured as a function of the transverse momentum and rapidity of both the hadronic and leptonic pseudo-top-quark as well as the transverse momentum, rapidity and invariant mass of the pseudo-top-quark pair system. The measurements are corrected for detector effects and are presented within a kinematic range that closely matches the detector acceptance.

Source code: ATLAS_2015_I1345452.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
 10namespace Rivet {
 11
 12
 13  /// @brief ATLAS 7 TeV pseudo-top analysis
 14  ///
 15  /// @author K .Finelli <kevin.finelli@cern.ch>
 16  /// @author A. Saavedra <a.saavedra@physics.usyd.edu.au>
 17  /// @author L. Lan <llan@physics.usyd.edu.au>
 18  class ATLAS_2015_I1345452 : public Analysis {
 19  public:
 20
 21    /// Constructor
 22    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2015_I1345452);
 23
 24
 25    void init() {
 26      // Eta ranges
 27      Cut eta_full = Cuts::abseta < 5.0 && Cuts::pT >= 1.0*MeV;
 28      Cut eta_lep = Cuts::abseta < 2.5;
 29
 30      // All final state particles
 31      FinalState fs(eta_full);
 32
 33      // Get photons to dress leptons
 34      IdentifiedFinalState photons(fs);
 35      photons.acceptIdPair(PID::PHOTON);
 36
 37      // Projection to find the electrons
 38      IdentifiedFinalState el_id(fs);
 39      el_id.acceptIdPair(PID::ELECTRON);
 40
 41      PromptFinalState electrons(el_id);
 42      electrons.acceptTauDecays(true);
 43      declare(electrons, "electrons");
 44
 45      LeptonFinder dressedelectrons(electrons, photons, 0.1, eta_lep && Cuts::pT > 25*GeV);
 46      declare(dressedelectrons, "dressedelectrons");
 47
 48      LeptonFinder ewdressedelectrons(electrons, photons, 0.1, eta_full);
 49      declare(ewdressedelectrons, "ewdressedelectrons");
 50
 51      LeptonFinder vetodressedelectrons(electrons, photons, 0.1, eta_lep && Cuts::pT > 15*GeV);
 52      declare(vetodressedelectrons, "vetodressedelectrons");
 53
 54      // Projection to find the muons
 55      IdentifiedFinalState mu_id(fs);
 56      mu_id.acceptIdPair(PID::MUON);
 57      PromptFinalState muons(mu_id);
 58      muons.acceptTauDecays(true);
 59      declare(muons, "muons");
 60      LeptonFinder dressedmuons(muons, photons, 0.1, eta_lep && Cuts::pT > 25*GeV);
 61      declare(dressedmuons, "dressedmuons");
 62      LeptonFinder ewdressedmuons(muons, photons, 0.1, eta_full);
 63      declare(ewdressedmuons, "ewdressedmuons");
 64      LeptonFinder vetodressedmuons(muons, photons, 0.1, eta_lep && Cuts::pT > 15*GeV);
 65      declare(vetodressedmuons, "vetodressedmuons");
 66
 67      // Projection to find neutrinos and produce MET
 68      IdentifiedFinalState nu_id;
 69      nu_id.acceptNeutrinos();
 70      PromptFinalState neutrinos(nu_id);
 71      neutrinos.acceptTauDecays(true);
 72      declare(neutrinos, "neutrinos");
 73
 74      // Jet clustering.
 75      VetoedFinalState vfs;
 76      vfs.addVetoOnThisFinalState(ewdressedelectrons);
 77      vfs.addVetoOnThisFinalState(ewdressedmuons);
 78      vfs.addVetoOnThisFinalState(neutrinos);
 79      FastJets jets(vfs, JetAlg::ANTIKT, 0.4);
 80      jets.useInvisibles();
 81      declare(jets, "jets");
 82
 83      //pseudotop leptons and hadrons
 84      book(_h["ptpseudotophadron_mu"]     , 1, 1, 2);
 85      book(_h["ptpseudotophadron_el"]     , 2, 1, 2);
 86      book(_h["absrappseudotophadron_mu"] , 3, 1, 2);
 87      book(_h["absrappseudotophadron_el"] , 4, 1, 2);
 88      book(_h["ptpseudotoplepton_mu"]     , 5, 1, 2);
 89      book(_h["ptpseudotoplepton_el"]     , 6, 1, 2);
 90      book(_h["absrappseudotoplepton_mu"] , 7, 1, 2);
 91      book(_h["absrappseudotoplepton_el"] , 8, 1, 2);
 92      book(_h["ptttbar_mu"]               , 9, 1, 2);
 93      book(_h["ptttbar_el"]               ,10, 1, 2);
 94      book(_h["absrapttbar_mu"]           ,11, 1, 2);
 95      book(_h["absrapttbar_el"]           ,12, 1, 2);
 96      book(_h["ttbarmass_mu"]             ,13, 1, 2);
 97      book(_h["ttbarmass_el"]             ,14, 1, 2);
 98      book(_h["ptpseudotophadron"]        ,15, 1, 2);
 99      book(_h["absrappseudotophadron"]    ,16, 1, 2);
100      book(_h["ptpseudotoplepton"]        ,17, 1, 2);
101      book(_h["absrappseudotoplepton"]    ,18, 1, 2);
102      book(_h["ptttbar"]                  ,19, 1, 2);
103      book(_h["absrapttbar"]              ,20, 1, 2);
104      book(_h["ttbarmass"]                ,21, 1, 2);
105
106    }
107
108    void analyze(const Event& event) {
109
110      // Get the selected objects, using the projections.
111      _dressedelectrons     = apply<LeptonFinder>(  event, "dressedelectrons").dressedLeptons();
112      _vetodressedelectrons = apply<LeptonFinder>(  event, "vetodressedelectrons").dressedLeptons();
113      _dressedmuons         = apply<LeptonFinder>(  event, "dressedmuons").dressedLeptons();
114      _vetodressedmuons     = apply<LeptonFinder>(  event, "vetodressedmuons").dressedLeptons();
115      _neutrinos            = apply<PromptFinalState>(event, "neutrinos").particlesByPt();
116      const Jets& all_jets  = apply<FastJets>(        event, "jets").jetsByPt(Cuts::pT > 25.0*GeV && Cuts::abseta < 2.5);
117
118      //get true l+jets events by removing events with more than 1 electron||muon neutrino
119      unsigned int n_elmu_neutrinos = 0;
120      for (const Particle & p : _neutrinos) {
121        if (p.abspid() == 12 || p.abspid() == 14)  ++n_elmu_neutrinos;
122      }
123      if (n_elmu_neutrinos != 1)  vetoEvent;
124
125      DressedLepton *lepton;
126      if ( _dressedelectrons.size())  lepton = &_dressedelectrons[0];
127      else if (_dressedmuons.size())  lepton = &_dressedmuons[0];
128      else vetoEvent;
129
130      // Calculate the missing ET, using the prompt neutrinos only (really?)
131      /// @todo Why not use MissingMomentum?
132      FourMomentum met;
133      for (const Particle& p : _neutrinos)  met += p.momentum();
134
135      //remove jets if they are within dR < 0.2 of lepton
136      Jets jets;
137      for(const Jet& jet : all_jets) {
138        bool keep = true;
139        for (const DressedLepton& el : _vetodressedelectrons) {
140          keep &= deltaR(jet, el) >= 0.2;
141        }
142        if (keep)  jets += jet;
143      }
144
145      bool overlap = false;
146      Jets bjets, lightjets;
147      for (unsigned int i = 0; i < jets.size(); ++i) {
148        const Jet& jet = jets[i];
149        for (const DressedLepton& el : _dressedelectrons)  overlap |= deltaR(jet, el) < 0.4;
150        for (const DressedLepton& mu : _dressedmuons)      overlap |= deltaR(jet, mu) < 0.4;
151        for (unsigned int j = i + 1; j < jets.size(); ++j) {
152          overlap |= deltaR(jet, jets[j]) < 0.5;
153        }
154        //// Count the number of b-tags
155        bool b_tagged = false;           //  This is closer to the
156        Particles bTags = jet.bTags();   //  analysis. Something
157        for ( Particle b : bTags ) {  //  about ghost-associated
158          b_tagged |= b.pT() > 5*GeV;    //  B-hadrons
159        }                                //
160        if ( b_tagged )  bjets += jet;
161        else lightjets += jet;
162      }
163
164      // remove events with object overlap
165      if (overlap) vetoEvent;
166
167      if (bjets.size() < 2 || lightjets.size() < 2)  vetoEvent;
168
169      FourMomentum pbjet1; //Momentum of bjet1
170      FourMomentum pbjet2; //Momentum of bjet2
171      if ( deltaR(bjets[0], *lepton) <= deltaR(bjets[1], *lepton) ) {
172        pbjet1 = bjets[0].momentum();
173        pbjet2 = bjets[1].momentum();
174      } else {
175        pbjet1 = bjets[1].momentum();
176        pbjet2 = bjets[0].momentum();
177      }
178
179      FourMomentum pjet1; // Momentum of jet1
180      if (lightjets.size())  pjet1 = lightjets[0].momentum();
181
182      FourMomentum pjet2; // Momentum of jet 2
183      if (lightjets.size() > 1)  pjet2 = lightjets[1].momentum();
184
185      double pz = computeneutrinoz(lepton->momentum(), met);
186      FourMomentum ppseudoneutrino( sqrt(sqr(met.px()) + sqr(met.py()) + sqr(pz)), met.px(), met.py(), pz);
187
188      //compute leptonic, hadronic, combined pseudo-top
189      FourMomentum ppseudotoplepton = lepton->momentum() + ppseudoneutrino + pbjet1;
190      FourMomentum ppseudotophadron = pbjet2 + pjet1 + pjet2;
191      FourMomentum pttbar = ppseudotoplepton + ppseudotophadron;
192
193      // Evaluate basic event selection
194      bool pass_eljets = (_dressedelectrons.size() == 1) &&
195                                                                                           (_vetodressedelectrons.size() < 2) &&
196        (_vetodressedmuons.empty()) &&
197        (met.pT() > 30*GeV) &&
198                                                                                           (_mT(_dressedelectrons[0].momentum(), met) > 35*GeV) &&
199                                                                                           (jets.size() >= 4);
200      bool pass_mujets = (_dressedmuons.size() == 1) &&
201        (_vetodressedmuons.size() < 2) &&
202        (_vetodressedelectrons.empty()) &&
203        (met.pT() > 30*GeV) &&
204        (_mT(_dressedmuons[0].momentum(), met) > 35*GeV) &&
205        (jets.size() >= 4);
206
207      // basic event selection requirements
208      if (!pass_eljets && !pass_mujets) vetoEvent;
209
210      // Fill histograms
211      //pseudotop hadrons and leptons fill histogram
212      _h["ptpseudotoplepton"]->fill(    ppseudotoplepton.pt()); //pT of pseudo top lepton
213      _h["absrappseudotoplepton"]->fill(ppseudotoplepton.absrap()); //absolute rapidity of pseudo top lepton
214      _h["ptpseudotophadron"]->fill(    ppseudotophadron.pt()); //pT of pseudo top hadron
215      _h["absrappseudotophadron"]->fill(ppseudotophadron.absrap()); //absolute rapidity of pseudo top hadron
216      _h["absrapttbar"]->fill(          pttbar.absrap()); //absolute rapidity of ttbar
217      _h["ttbarmass"]->fill(            pttbar.mass()); //mass of ttbar
218      _h["ptttbar"]->fill(              pttbar.pt()); //fill pT of ttbar in combined channel
219
220      if (pass_eljets) { // electron channel fill histogram
221        _h["ptpseudotoplepton_el"]->fill(    ppseudotoplepton.pt()); //pT of pseudo top lepton
222        _h["absrappseudotoplepton_el"]->fill(ppseudotoplepton.absrap()); //absolute rapidity of pseudo top lepton
223        _h["ptpseudotophadron_el"]->fill(    ppseudotophadron.pt()); //pT of pseudo top hadron
224        _h["absrappseudotophadron_el"]->fill(ppseudotophadron.absrap()); //absolute rapidity of pseudo top hadron
225        _h["absrapttbar_el"]->fill(          pttbar.absrap()); //absolute rapidity of ttbar
226        _h["ttbarmass_el"]->fill(            pttbar.mass()); // fill electron channel ttbar mass
227        _h["ptttbar_el"]->fill(              pttbar.pt()); //fill pT of ttbar in electron channel
228      }
229      else { // muon channel fill histogram
230        _h["ptpseudotoplepton_mu"]->fill(    ppseudotoplepton.pt()); //pT of pseudo top lepton
231        _h["absrappseudotoplepton_mu"]->fill(ppseudotoplepton.absrap()); //absolute rapidity of pseudo top lepton
232        _h["ptpseudotophadron_mu"]->fill(    ppseudotophadron.pt()); //pT of pseudo top hadron
233        _h["absrappseudotophadron_mu"]->fill(ppseudotophadron.absrap()); //absolute rapidity of pseudo top hadron
234        _h["absrapttbar_mu"]->fill(          pttbar.absrap()); //absolute rapidity of ttbar
235        _h["ttbarmass_mu"]->fill(            pttbar.mass()); //fill muon channel histograms
236        _h["ptttbar_mu"]->fill(              pttbar.pt()); //fill pT of ttbar in electron channel
237      }
238    }
239
240    void finalize() {
241      // Normalize to cross-section
242      const double scalefactor(crossSection()/picobarn / sumOfWeights());
243      for (map<string, Histo1DPtr>::iterator hit = _h.begin(); hit != _h.end(); ++hit) {
244        double sf = scalefactor;
245        if ( (hit->first).find("_") == std::string::npos )  sf *= 0.5;
246        scale(hit->second, sf);
247      }
248    }
249
250  private:
251
252
253    double computeneutrinoz(const FourMomentum& lepton, FourMomentum& met) const {
254      //computing z component of neutrino momentum given lepton and met
255      double pzneutrino;
256      double m_W = 80.399; // in GeV, given in the paper
257      double k = (( sqr( m_W ) - sqr( lepton.mass() ) ) / 2 ) + (lepton.px() * met.px() + lepton.py() * met.py());
258      double a = sqr ( lepton.E() )- sqr ( lepton.pz() );
259      double b = -2*k*lepton.pz();
260      double c = sqr( lepton.E() ) * sqr( met.pT() ) - sqr( k );
261      double discriminant = sqr(b) - 4 * a * c;
262      double quad[2] = { (- b - sqrt(discriminant)) / (2 * a), (- b + sqrt(discriminant)) / (2 * a) }; //two possible quadratic solns
263      if (discriminant < 0)  pzneutrino = - b / (2 * a); //if the discriminant is negative
264      else { //if the discriminant is greater than or equal to zero, take the soln with smallest absolute value
265        double absquad[2];
266        for (int n=0; n<2; ++n)  absquad[n] = fabs(quad[n]);
267        if (absquad[0] < absquad[1])  pzneutrino = quad[0];
268        else                          pzneutrino = quad[1];
269      }
270      if ( !std::isfinite(pzneutrino) )  std::cout << "Found non-finite value\n";
271      return pzneutrino;
272    }
273
274    double _mT(const FourMomentum &l, FourMomentum &nu) const {
275      return sqrt( 2 * l.pT() * nu.pT() * (1 - cos(deltaPhi(l, nu))) );
276    }
277
278    /// @name Objects that are used by the event selection decisions
279    DressedLeptons _dressedelectrons, _vetodressedelectrons, _dressedmuons, _vetodressedmuons;
280    Particles _neutrinos;
281    map<string, Histo1DPtr> _h;
282  };
283
284  RIVET_DECLARE_PLUGIN(ATLAS_2015_I1345452);
285
286}