rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

D0_1996_S3214044

Topological distributions of inclusive three- and four-jet events
Experiment: D0 (Tevatron Run 1)
Inspire ID: 399364
Status: VALIDATED
Authors:
  • Frank Siegert
References: Beams: p- p+
Beam energies: (900.0, 900.0) GeV
Run details:
  • $p \bar{p} \to$ jets at 1800 GeV with minimum jet pT in analysis = 20 GeV

The global topologies of inclusive three- and four-jet events produced in pbar p interactions are described. The three- and four-jet events are selected from data recorded by the D0 detector at the Fermilab Tevatron Collider operating at a center-of-mass energy of $\sqrt{s}$=1800 GeV. The studies also show that the topological distributions of the different subprocesses involving different numbers of quarks are very similar and reproduce the measured distributions well. The parton-shower Monte Carlo generators provide a less satisfactory description of the topologies of the three- and four-jet events.

Source code: D0_1996_S3214044.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Math/LorentzTrans.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// D0 topological distributions of 3- and 4-jet events.
 11  class D0_1996_S3214044 : public Analysis {
 12  public:
 13
 14    RIVET_DEFAULT_ANALYSIS_CTOR(D0_1996_S3214044);
 15
 16
 17    /// @name Analysis methods
 18    /// @{
 19
 20    /// Book histograms
 21    void init() {
 22      const FinalState fs;
 23      declare(fs, "FS");
 24      /// @todo Use correct jet algorithm --- tried FJ3 D0RunICone but does
 25      // not look as good as the Run2 cone alg used here
 26      declare(FastJets(fs, FastJets::D0ILCONE, 0.7), "ConeJets");
 27
 28      book(_h_3j_x3 ,1, 1, 1);
 29      book(_h_3j_x5 ,2, 1, 1);
 30      book(_h_3j_costheta3 ,3, 1, 1);
 31      book(_h_3j_psi ,4, 1, 1);
 32      book(_h_3j_mu34 ,5, 1, 1);
 33      book(_h_3j_mu35 ,6, 1, 1);
 34      book(_h_3j_mu45 ,7, 1, 1);
 35
 36      book(_h_4j_x3 ,8, 1, 1);
 37      book(_h_4j_x4 ,9, 1, 1);
 38      book(_h_4j_x5 ,10, 1, 1);
 39      book(_h_4j_x6 ,11, 1, 1);
 40      book(_h_4j_costheta3 ,12, 1, 1);
 41      book(_h_4j_costheta4 ,13, 1, 1);
 42      book(_h_4j_costheta5 ,14, 1, 1);
 43      book(_h_4j_costheta6 ,15, 1, 1);
 44      book(_h_4j_cosomega34 ,16, 1, 1);
 45      book(_h_4j_cosomega35 ,17, 1, 1);
 46      book(_h_4j_cosomega36 ,18, 1, 1);
 47      book(_h_4j_cosomega45 ,19, 1, 1);
 48      book(_h_4j_cosomega46 ,20, 1, 1);
 49      book(_h_4j_cosomega56 ,21, 1, 1);
 50      book(_h_4j_mu34 ,22, 1, 1);
 51      book(_h_4j_mu35 ,23, 1, 1);
 52      book(_h_4j_mu36 ,24, 1, 1);
 53      book(_h_4j_mu45 ,25, 1, 1);
 54      book(_h_4j_mu46 ,26, 1, 1);
 55      book(_h_4j_mu56 ,27, 1, 1);
 56      book(_h_4j_theta_BZ ,28, 1, 1);
 57      book(_h_4j_costheta_NR ,29, 1, 1);
 58
 59    }
 60
 61
 62    void analyze(const Event& event) {
 63      Jets jets_in = apply<FastJets>(event, "ConeJets")
 64        .jets(Cuts::Et > 20*GeV && Cuts::abseta < 3, cmpMomByEt);
 65
 66      Jets jets_isolated;
 67      for (size_t i = 0; i < jets_in.size(); ++i) {
 68        bool isolated = true;
 69        for (size_t j = 0; j < jets_in.size(); ++j) {
 70          if (i != j && deltaR(jets_in[i], jets_in[j]) < 1.4) {
 71            isolated = false;
 72            break;
 73          }
 74        }
 75        if (isolated) jets_isolated.push_back(jets_in[i]);
 76      }
 77
 78      if (jets_isolated.size() == 0 || jets_isolated[0].Et() < 60.0*GeV) vetoEvent;
 79
 80      if (jets_isolated.size() > 2) _threeJetAnalysis(jets_isolated);
 81      if (jets_isolated.size() > 3) _fourJetAnalysis(jets_isolated);
 82    }
 83
 84
 85    void finalize() {
 86      normalize(_h_3j_x3, 1.0);
 87      normalize(_h_3j_x5, 1.0);
 88      normalize(_h_3j_costheta3, 1.0);
 89      normalize(_h_3j_psi, 1.0);
 90      normalize(_h_3j_mu34, 1.0);
 91      normalize(_h_3j_mu35, 1.0);
 92      normalize(_h_3j_mu45, 1.0);
 93      normalize(_h_4j_x3, 1.0);
 94      normalize(_h_4j_x4, 1.0);
 95      normalize(_h_4j_x5, 1.0);
 96      normalize(_h_4j_x6, 1.0);
 97      normalize(_h_4j_costheta3, 1.0);
 98      normalize(_h_4j_costheta4, 1.0);
 99      normalize(_h_4j_costheta5, 1.0);
100      normalize(_h_4j_costheta6, 1.0);
101      normalize(_h_4j_cosomega34, 1.0);
102      normalize(_h_4j_cosomega35, 1.0);
103      normalize(_h_4j_cosomega36, 1.0);
104      normalize(_h_4j_cosomega45, 1.0);
105      normalize(_h_4j_cosomega46, 1.0);
106      normalize(_h_4j_cosomega56, 1.0);
107      normalize(_h_4j_mu34, 1.0);
108      normalize(_h_4j_mu35, 1.0);
109      normalize(_h_4j_mu36, 1.0);
110      normalize(_h_4j_mu45, 1.0);
111      normalize(_h_4j_mu46, 1.0);
112      normalize(_h_4j_mu56, 1.0);
113      normalize(_h_4j_theta_BZ, 1.0);
114      normalize(_h_4j_costheta_NR, 1.0);
115    }
116
117    /// @}
118
119
120  private:
121
122    /// @name Helper functions
123    /// @{
124
125    void _threeJetAnalysis(const Jets& jets) {
126      // >=3 jet events
127      FourMomentum jjj(jets[0].momentum()+jets[1].momentum()+jets[2].momentum());
128      const double sqrts = _safeMass(jjj);
129      if (sqrts<200*GeV) {
130        return;
131      }
132
133      const LorentzTransform cms_boost = LorentzTransform::mkFrameTransformFromBeta(jjj.betaVec());
134      vector<FourMomentum> jets_boosted;
135      for (Jet jet : jets) {
136        jets_boosted.push_back(cms_boost.transform(jet.momentum()));
137      }
138      std::sort(jets_boosted.begin(), jets_boosted.end(), FourMomentum::byEDescending());
139      FourMomentum p3(jets_boosted[0]);
140      FourMomentum p4(jets_boosted[1]);
141      FourMomentum p5(jets_boosted[2]);
142
143      Vector3 beam1(0.0, 0.0, 1.0);
144      Vector3 p1xp3 = beam1.cross(p3.p3());
145      Vector3 p4xp5 = p4.p3().cross(p5.p3());
146      const double cospsi = p1xp3.dot(p4xp5)/p1xp3.mod()/p4xp5.mod();
147
148      _h_3j_x3->fill(2.0*p3.E()/sqrts);
149      _h_3j_x5->fill(2.0*p5.E()/sqrts);
150      _h_3j_costheta3->fill(fabs(cos(p3.theta())));
151      _h_3j_psi->fill(acos(cospsi)/degree);
152      _h_3j_mu34->fill(_safeMass(FourMomentum(p3+p4))/sqrts);
153      _h_3j_mu35->fill(_safeMass(FourMomentum(p3+p5))/sqrts);
154      _h_3j_mu45->fill(_safeMass(FourMomentum(p4+p5))/sqrts);
155    }
156
157
158    void _fourJetAnalysis(const Jets& jets) {
159      // >=4 jet events
160      FourMomentum jjjj(jets[0].momentum() + jets[1].momentum() + jets[2].momentum()+ jets[3].momentum());
161      const double sqrts = _safeMass(jjjj);
162      if (sqrts < 200*GeV) return;
163
164      const LorentzTransform cms_boost = LorentzTransform::mkFrameTransformFromBeta(jjjj.betaVec());
165      vector<FourMomentum> jets_boosted;
166      for (Jet jet : jets) {
167        jets_boosted.push_back(cms_boost.transform(jet.momentum()));
168      }
169      sort(jets_boosted.begin(), jets_boosted.end(), FourMomentum::byEDescending());
170      FourMomentum p3(jets_boosted[0]);
171      FourMomentum p4(jets_boosted[1]);
172      FourMomentum p5(jets_boosted[2]);
173      FourMomentum p6(jets_boosted[3]);
174
175      Vector3 p3xp4 = p3.p3().cross(p4.p3());
176      Vector3 p5xp6 = p5.p3().cross(p6.p3());
177      const double costheta_BZ = p3xp4.dot(p5xp6)/p3xp4.mod()/p5xp6.mod();
178      const double costheta_NR = (p3.p3()-p4.p3()).dot(p5.p3()-p6.p3())/
179        (p3.p3()-p4.p3()).mod()/(p5.p3()-p6.p3()).mod();
180
181      _h_4j_x3->fill(2.0*p3.E()/sqrts);
182      _h_4j_x4->fill(2.0*p4.E()/sqrts);
183      _h_4j_x5->fill(2.0*p5.E()/sqrts);
184      _h_4j_x6->fill(2.0*p6.E()/sqrts);
185      _h_4j_costheta3->fill(fabs(cos(p3.theta())));
186      _h_4j_costheta4->fill(fabs(cos(p4.theta())));
187      _h_4j_costheta5->fill(fabs(cos(p5.theta())));
188      _h_4j_costheta6->fill(fabs(cos(p6.theta())));
189      _h_4j_cosomega34->fill(cos(p3.angle(p4)));
190      _h_4j_cosomega35->fill(cos(p3.angle(p5)));
191      _h_4j_cosomega36->fill(cos(p3.angle(p6)));
192      _h_4j_cosomega45->fill(cos(p4.angle(p5)));
193      _h_4j_cosomega46->fill(cos(p4.angle(p6)));
194      _h_4j_cosomega56->fill(cos(p5.angle(p6)));
195      _h_4j_mu34->fill(_safeMass(FourMomentum(p3+p4))/sqrts);
196      _h_4j_mu35->fill(_safeMass(FourMomentum(p3+p5))/sqrts);
197      _h_4j_mu36->fill(_safeMass(FourMomentum(p3+p6))/sqrts);
198      _h_4j_mu45->fill(_safeMass(FourMomentum(p4+p5))/sqrts);
199      _h_4j_mu46->fill(_safeMass(FourMomentum(p4+p6))/sqrts);
200      _h_4j_mu56->fill(_safeMass(FourMomentum(p5+p6))/sqrts);
201      _h_4j_theta_BZ->fill(acos(fabs(costheta_BZ))/degree);
202      _h_4j_costheta_NR->fill(fabs(costheta_NR));
203
204    }
205
206    double _safeMass(const FourMomentum& p) {
207      double mass2=p.mass2();
208      if (mass2>0.0) return sqrt(mass2);
209      else if (mass2<-1.0e-5) {
210        MSG_WARNING("m2 = " << m2 << ". Assuming m2=0.");
211        return 0.0;
212      }
213      else return 0.0;
214    }
215
216    /// @}
217
218
219  private:
220
221    /// @name Histograms
222    /// @{
223    Histo1DPtr _h_3j_x3;
224    Histo1DPtr _h_3j_x5;
225    Histo1DPtr _h_3j_costheta3;
226    Histo1DPtr _h_3j_psi;
227    Histo1DPtr _h_3j_mu34;
228    Histo1DPtr _h_3j_mu35;
229    Histo1DPtr _h_3j_mu45;
230
231    Histo1DPtr _h_4j_x3;
232    Histo1DPtr _h_4j_x4;
233    Histo1DPtr _h_4j_x5;
234    Histo1DPtr _h_4j_x6;
235    Histo1DPtr _h_4j_costheta3;
236    Histo1DPtr _h_4j_costheta4;
237    Histo1DPtr _h_4j_costheta5;
238    Histo1DPtr _h_4j_costheta6;
239    Histo1DPtr _h_4j_cosomega34;
240    Histo1DPtr _h_4j_cosomega35;
241    Histo1DPtr _h_4j_cosomega36;
242    Histo1DPtr _h_4j_cosomega45;
243    Histo1DPtr _h_4j_cosomega46;
244    Histo1DPtr _h_4j_cosomega56;
245    Histo1DPtr _h_4j_mu34;
246    Histo1DPtr _h_4j_mu35;
247    Histo1DPtr _h_4j_mu36;
248    Histo1DPtr _h_4j_mu45;
249    Histo1DPtr _h_4j_mu46;
250    Histo1DPtr _h_4j_mu56;
251    Histo1DPtr _h_4j_theta_BZ;
252    Histo1DPtr _h_4j_costheta_NR;
253    /// @}
254
255  };
256
257
258
259  RIVET_DECLARE_ALIASED_PLUGIN(D0_1996_S3214044, D0_1996_I399364);
260
261}