rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1310835

H(125) -> 4l at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1310835
Status: VALIDATED
Authors:
  • Jonathan Stahlman
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • p + p -> H (-> 4 l) + X at 8 TeV

Measurements of fiducial and differential cross sections of Higgs boson production in the $H\to ZZ^\ast\to 4\ell$ decay channel are presented. The cross sections are determined within a fiducial phase space and corrected for detection efficiency and resolution effects. They are based on 20.3 fb$^{-1}$ of $pp$ collision data, produced at $\sqrt{s}=8$ TeV centre-of-mass energy at the LHC and recorded by the ATLAS detector. The differential measurements are performed in bins of transverse momentum and rapidity of the four-lepton system, the invariant mass of the subleading lepton pair and the decay angle of the leading lepton pair with respect to the beam line in the four-lepton rest frame, as well as the number of jets and the transverse momentum of the leading jet. The measured cross sections are compared to selected theoretical calculations of the Standard Model expectations. No significant deviation from any of the tested predictions is found.

Source code: ATLAS_2014_I1310835.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6#include "Rivet/Projections/DressedLeptons.hh"
  7
  8namespace Rivet {
  9
 10  /// @brief H(125)->ZZ->4l at 8 TeV
 11  class ATLAS_2014_I1310835 : public Analysis {
 12  public:
 13
 14    /// Default constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1310835);
 16
 17    void init() {
 18      const FinalState fs(Cuts::abseta < 5.0);
 19
 20      PromptFinalState photons(Cuts::abspid == PID::PHOTON);
 21
 22      PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON);
 23
 24      PromptFinalState bare_mu(Cuts::abspid == PID::MUON);
 25
 26      // Selection: lepton selection
 27      Cut etaranges_el = Cuts::abseta < 2.47 && Cuts::pT > 7*GeV; 
 28      DressedLeptons electron_sel4l(photons, bare_el, 0.1, etaranges_el, false);
 29      declare(electron_sel4l, "electrons");
 30 
 31      Cut etaranges_mu = Cuts::abseta < 2.7 && Cuts::pT > 6*GeV;
 32      DressedLeptons muon_sel4l(photons, bare_mu, 0.1, etaranges_mu, false);
 33      declare(muon_sel4l, "muons");
 34
 35      FastJets jetpro(fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::NONE, JetAlg::Invisibles::NONE);
 36      declare(jetpro, "jet");
 37
 38      // Book histos
 39      book(_h_pt          , 1, 1, 1);
 40      book(_h_rapidity    , 2, 1, 1);
 41      book(_h_m34         , 3, 1, 1);
 42      book(_h_costheta    , 4, 1, 1);
 43      book(_h_njets       , 5, 1, 1);
 44      book(_h_leadingjetpt, 6, 1, 1);
 45
 46    }
 47
 48
 49
 50    /// Do the analysis
 51    void analyze(const Event& e) {
 52      
 53      ////////////////////////////////////////////////////////////////////
 54      // preselection of leptons for ZZ-> llll final state
 55      ////////////////////////////////////////////////////////////////////
 56
 57      const vector<DressedLepton>& mu_sel4l = applyProjection<DressedLeptons>(e, "muons").dressedLeptons();
 58      const vector<DressedLepton>& el_sel4l = applyProjection<DressedLeptons>(e, "electrons").dressedLeptons();
 59
 60      vector<DressedLepton> leptonsFS_sel4l;
 61      leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), mu_sel4l.begin(), mu_sel4l.end() );
 62      leptonsFS_sel4l.insert( leptonsFS_sel4l.end(), el_sel4l.begin(), el_sel4l.end() );
 63
 64      /////////////////////////////////////////////////////////////////////////////
 65      /// H->ZZ->4l pairing
 66      /////////////////////////////////////////////////////////////////////////////
 67 
 68      size_t el_p = 0;
 69      size_t el_n = 0;
 70      size_t mu_p = 0; 
 71      size_t mu_n = 0;
 72      
 73      for (const Particle& l : leptonsFS_sel4l) {
 74        if (l.abspid() == PID::ELECTRON) {
 75          if (l.pid() < 0)  ++el_n;
 76          if (l.pid() > 0)  ++el_p;
 77        }
 78        else if (l.abspid() == PID::MUON) {
 79          if (l.pid() < 0)  ++mu_n;
 80          if (l.pid() > 0)  ++mu_p;
 81        }
 82      }
 83            
 84      bool pass_sfos = ( (el_p >=2 && el_n >=2) || (mu_p >=2 && mu_n >=2) || (el_p >=1 && el_n >=1 && mu_p >=1 && mu_n >=1) );
 85      
 86      if (!pass_sfos)  vetoEvent;
 87
 88      Zstate Z1, Z2, Zcand;
 89      size_t n_parts = leptonsFS_sel4l.size();
 90      size_t l1_index = 0;
 91      size_t l2_index = 0;
 92
 93      // determine Z1 first
 94      double min_mass_diff = -1;
 95      for (size_t i = 0; i < n_parts; ++i) {
 96        for (size_t j = 0; j < n_parts; ++j) {
 97          if (i >= j)  continue;
 98
 99          if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid())  continue; //only pair SFOS leptons
100
101          Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
102          double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
103         
104          if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
105            min_mass_diff = mass_diff;
106            Z1 = Zcand;
107            l1_index = i;
108            l2_index = j;
109          }
110        }
111      }
112
113      //determine Z2 second
114      min_mass_diff = -1;
115      for (size_t i = 0; i < n_parts; ++i) {
116        if (i == l1_index || i == l2_index)  continue;
117        for (size_t j = 0; j < n_parts; ++j) {
118          if (j == l1_index || j == l2_index || i >= j)  continue;
119
120          if (leptonsFS_sel4l[i].pid() != -1*leptonsFS_sel4l[j].pid())  continue; // only pair SFOS leptons
121
122          Zcand = Zstate( ParticlePair(leptonsFS_sel4l[i], leptonsFS_sel4l[j]) );
123          double mass_diff = fabs( Zcand.mom().mass() - 91.1876 );
124
125          if (min_mass_diff == -1 || mass_diff < min_mass_diff) {
126            min_mass_diff = mass_diff;
127            Z2 = Zcand;
128          }
129        }
130      }
131
132      Particles leptons_sel4l;
133      leptons_sel4l.push_back(Z1.first);
134      leptons_sel4l.push_back(Z1.second);
135      leptons_sel4l.push_back(Z2.first);
136      leptons_sel4l.push_back(Z2.second);
137
138      ////////////////////////////////////////////////////////////////////////////
139      // Kinematic Requirements
140      ///////////////////////////////////////////////////////////////////////////
141      
142      //leading lepton pT requirement
143      std::vector<double> lepton_pt;
144      for (const Particle& i : leptons_sel4l) lepton_pt.push_back(i.pT() / GeV);
145      std::sort(lepton_pt.begin(), lepton_pt.end(), [](const double pT1, const double pT2) { return pT1 > pT2; });
146      
147      if (!(lepton_pt[0] > 20*GeV && lepton_pt[1] > 15*GeV && lepton_pt[2] > 10*GeV))  vetoEvent;
148      
149      //invariant mass requirements
150      if (!(inRange(Z1.mom().mass(), 50*GeV, 106*GeV) && inRange(Z2.mom().mass(), 12*GeV, 115*GeV)))  vetoEvent;
151      
152      //lepton separation requirements
153      for (unsigned int i = 0; i < 4; ++i) {
154        for (unsigned int j = 0; j < 4; ++j) {
155          if (i >= j) continue;
156          double dR = deltaR(leptons_sel4l[i], leptons_sel4l[j]);
157          bool sameflavor = leptons_sel4l[i].abspid() == leptons_sel4l[j].abspid();
158
159          if ( sameflavor && dR < 0.1)  vetoEvent;
160          if (!sameflavor && dR < 0.2)  vetoEvent;
161        }
162      }
163
164      // J/Psi veto requirement
165      for (unsigned int i = 0; i < 4; ++i) {
166        for (unsigned int j = 0; j < 4; ++j) {
167          if (i >= j) continue;
168          if ( leptons_sel4l[i].pid() != -1*leptons_sel4l[j].pid() )  continue;
169          if ((leptons_sel4l[i].momentum() + leptons_sel4l[j].momentum()).mass() <= 5*GeV)  vetoEvent;
170        }
171      }
172 
173      // 4-lepton invariant mass requirement
174      double m4l = (Z1.mom() + Z2.mom()).mass();
175      if (!(inRange(m4l, 118*GeV, 129*GeV)))  vetoEvent;
176  
177  
178      ////////////////////////////////////////////////////////////////////////////
179      // Higgs observables
180      ///////////////////////////////////////////////////////////////////////////
181      FourMomentum Higgs = Z1.mom() + Z2.mom();
182
183      double H4l_pt       = Higgs.pt()/GeV; 
184      double H4l_rapidity = Higgs.absrap(); 
185      LorentzTransform HRF_boost;
186      //HRF_boost.mkFrameTransformFromBeta(Higgs.betaVec());
187      HRF_boost.setBetaVec(- Higgs.betaVec());
188      FourMomentum Z1_in_HRF = HRF_boost.transform( Z1.mom() );
189      double H4l_costheta = fabs(cos( Z1_in_HRF.theta())); 
190      double H4l_m34      = Z2.mom().mass()/GeV;
191      
192      ////////////////////////////////////////////////////////////////////////////
193      // Jet observables
194      ///////////////////////////////////////////////////////////////////////////
195
196      Jets jets;
197      for (const Jet& jet : applyProjection<FastJets>(e, "jet").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.4)) {
198        bool overlaps = false;
199        for (const Particle& lep : leptonsFS_sel4l) {
200          if (lep.abspid() != PID::ELECTRON)  continue;
201          const double dR = deltaR(lep, jet);
202          if (dR < 0.2) { overlaps = true; break; }
203        }
204        if (!overlaps) jets += jet;
205      }
206      size_t n_jets = jets.size();
207      if (n_jets > 3)  n_jets = 3;
208
209      std::vector<double> jet_pt;
210      for (const Jet& i : jets) jet_pt.push_back(i.pT()/GeV);
211
212      double leading_jet_pt = n_jets? jet_pt[0] : 0.;
213
214      ////////////////////////////////////////////////////////////////////////////
215      // End of H->ZZ->llll selection: now fill histograms
216      ////////////////////////////////////////////////////////////////////////////
217
218
219      _h_pt->fill(H4l_pt);
220      _h_rapidity->fill(H4l_rapidity);
221      _h_costheta->fill(H4l_costheta);
222      _h_m34->fill(H4l_m34);
223      _h_njets->fill(n_jets + 1);
224      _h_leadingjetpt->fill(leading_jet_pt);
225
226
227    }
228
229
230    /// Generic Z candidate
231    struct Zstate : public ParticlePair {
232      Zstate() { }
233      Zstate(ParticlePair _particlepair) : ParticlePair(_particlepair) { }
234      FourMomentum mom() const { return first.momentum() + second.momentum(); }
235      operator FourMomentum() const { return mom(); }
236    };
237
238    /// Finalize
239    void finalize() {
240
241      const double norm = crossSection()/sumOfWeights()/femtobarn;
242      scale(_h_pt, norm);
243      scale(_h_rapidity, norm);
244      scale(_h_costheta, norm);
245      scale(_h_m34, norm);
246      scale(_h_njets, norm);
247      scale(_h_leadingjetpt, norm);
248    }
249
250
251  private:
252
253    Histo1DPtr _h_pt, _h_rapidity, _h_costheta;
254    Histo1DPtr _h_m34, _h_njets, _h_leadingjetpt;
255
256  };
257
258
259  // The hook for the plugin system
260  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1310835);
261
262}