rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

OPAL_2004_S6132243

Event shape distributions and moments in $e^+ e^-$ -> hadrons at 91--209 GeV
Experiment: OPAL (LEP 1 & 2)
Inspire ID: 669402
Status: VALIDATED
Authors:
  • Andy Buckley
References: Beams: e+ e-
Beam energies: (45.6, 45.6); (66.5, 66.5); (88.5, 88.5); (98.5, 98.5) GeV
Run details:
  • Hadronic $e^+ e^-$ events at 4 representative energies (91, 133, 177, 197). Runs need to have ISR suppressed, since the analysis was done using a cut of $\sqrt{s} - \sqrt{s_\text{reco}} < 1\,\text{GeV}$. Particles with a lifetime $> 3 \cdot 10^{-10}\,\text{s}$ are considered to be stable. Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples.

Measurement of $e^+ e^-$ event shape variable distributions and their 1st to 5th moments in LEP running from the Z pole to the highest LEP 2 energy of 209 GeV.Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples.

Source code: OPAL_2004_S6132243.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/ChargedFinalState.hh"
  6#include "Rivet/Projections/Sphericity.hh"
  7#include "Rivet/Projections/Thrust.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9#include "Rivet/Projections/ParisiTensor.hh"
 10#include "Rivet/Projections/Hemispheres.hh"
 11#include <cmath>
 12
 13namespace Rivet {
 14
 15
 16  /// @brief OPAL event shapes and moments at 91, 133, 177, and 197 GeV
 17  ///
 18  /// @author Andy Buckley
 19  class OPAL_2004_S6132243 : public Analysis {
 20  public:
 21
 22    RIVET_DEFAULT_ANALYSIS_CTOR(OPAL_2004_S6132243);
 23
 24
 25    /// @name Analysis methods
 26    /// @{
 27
 28    /// Energies: 91, 133, 177 (161-183), 197 (189-209) => index 0..4
 29    int getHistIndex() {
 30      int ih = -1;
 31      if (inRange(sqrtS()/GeV, 89.9, 91.5)) {
 32        ih = 0;
 33      } else if (isCompatibleWithSqrtS(133.)) {
 34        ih = 1;
 35      } else if (isCompatibleWithSqrtS(177.)) { // (161-183)
 36        ih = 2;
 37      } else if (isCompatibleWithSqrtS(197.)) { // (189-209)
 38        ih = 3;
 39      } else {
 40        stringstream ss;
 41        ss << "Invalid energy for OPAL_2004 analysis: "
 42           << sqrtS() << " GeV != 91, 133, 177, or 197 GeV";
 43        throw Error(ss.str());
 44      }
 45      assert(ih >= 0);
 46      return ih;
 47    }
 48
 49
 50    void init() {
 51      // Projections
 52      declare(Beam(), "Beams");
 53      const FinalState fs;
 54      declare(fs, "FS");
 55      const ChargedFinalState cfs;
 56      declare(cfs, "CFS");
 57      declare(FastJets(fs, FastJets::DURHAM, 0.7), "DurhamJets");
 58      declare(Sphericity(fs), "Sphericity");
 59      declare(ParisiTensor(fs), "Parisi");
 60      const Thrust thrust(fs);
 61      declare(thrust, "Thrust");
 62      declare(Hemispheres(thrust), "Hemispheres");
 63      _isqrts = getHistIndex();
 64
 65      // Book histograms
 66      book(_hist1MinusT[_isqrts]    ,1, 1, _isqrts+1);
 67      book(_histHemiMassH[_isqrts]  ,2, 1, _isqrts+1);
 68      book(_histCParam[_isqrts]     ,3, 1, _isqrts+1);
 69      book(_histHemiBroadT[_isqrts] ,4, 1, _isqrts+1);
 70      book(_histHemiBroadW[_isqrts] ,5, 1, _isqrts+1);
 71      book(_histY23Durham[_isqrts]  ,6, 1, _isqrts+1);
 72      book(_histTMajor[_isqrts]     ,7, 1, _isqrts+1);
 73      book(_histTMinor[_isqrts]     ,8, 1, _isqrts+1);
 74      book(_histAplanarity[_isqrts] ,9, 1, _isqrts+1);
 75      book(_histSphericity[_isqrts] ,10, 1, _isqrts+1);
 76      book(_histOblateness[_isqrts] ,11, 1, _isqrts+1);
 77      book(_histHemiMassL[_isqrts]  ,12, 1, _isqrts+1);
 78      book(_histHemiBroadN[_isqrts] ,13, 1, _isqrts+1);
 79      book(_histDParam[_isqrts]     ,14, 1, _isqrts+1);
 80      //
 81      book(_hist1MinusTMom[_isqrts]    ,15, 1, _isqrts+1);
 82      book(_histHemiMassHMom[_isqrts]  ,16, 1, _isqrts+1);
 83      book(_histCParamMom[_isqrts]     ,17, 1, _isqrts+1);
 84      book(_histHemiBroadTMom[_isqrts] ,18, 1, _isqrts+1);
 85      book(_histHemiBroadWMom[_isqrts] ,19, 1, _isqrts+1);
 86      book(_histY23DurhamMom[_isqrts]  ,20, 1, _isqrts+1);
 87      book(_histTMajorMom[_isqrts]     ,21, 1, _isqrts+1);
 88      book(_histTMinorMom[_isqrts]     ,22, 1, _isqrts+1);
 89      book(_histSphericityMom[_isqrts] ,23, 1, _isqrts+1);
 90      book(_histOblatenessMom[_isqrts] ,24, 1, _isqrts+1);
 91      book(_histHemiMassLMom[_isqrts]  ,25, 1, _isqrts+1);
 92      book(_histHemiBroadNMom[_isqrts] ,26, 1, _isqrts+1);
 93
 94      book(_sumWTrack2, "_sumWTrack2");
 95      book(_sumWJet3, "_sumWJet3");
 96    }
 97
 98
 99    void analyze(const Event& event) {
100      // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
101      const FinalState& cfs = apply<FinalState>(event, "CFS");
102      if (cfs.size() < 2) vetoEvent;
103
104      _sumWTrack2->fill();
105
106      // Thrusts
107      const Thrust& thrust = apply<Thrust>(event, "Thrust");
108      _hist1MinusT[_isqrts]->fill(1-thrust.thrust());
109      _histTMajor[_isqrts]->fill(thrust.thrustMajor());
110      _histTMinor[_isqrts]->fill(thrust.thrustMinor());
111      _histOblateness[_isqrts]->fill(thrust.oblateness());
112      for (int n = 1; n <= 5; ++n) {
113        _hist1MinusTMom[_isqrts]->fill(n, pow(1-thrust.thrust(), n));
114        _histTMajorMom[_isqrts]->fill(n, pow(thrust.thrustMajor(), n));
115        _histTMinorMom[_isqrts]->fill(n, pow(thrust.thrustMinor(), n));
116        _histOblatenessMom[_isqrts]->fill(n, pow(thrust.oblateness(), n));
117      }
118
119      // Jets
120      const FastJets& durjet = apply<FastJets>(event, "DurhamJets");
121      if (durjet.clusterSeq()) {
122        _sumWJet3->fill();
123        const double y23 = durjet.clusterSeq()->exclusive_ymerge_max(2);
124        if (y23>0.0) {
125          _histY23Durham[_isqrts]->fill(y23);
126          for (int n = 1; n <= 5; ++n) {
127            _histY23DurhamMom[_isqrts]->fill(n, pow(y23, n));
128          }
129        }
130      }
131
132      // Sphericities
133      const Sphericity& sphericity = apply<Sphericity>(event, "Sphericity");
134      const double sph = sphericity.sphericity();
135      const double apl = sphericity.aplanarity();
136      _histSphericity[_isqrts]->fill(sph);
137      _histAplanarity[_isqrts]->fill(apl);
138      for (int n = 1; n <= 5; ++n) {
139        _histSphericityMom[_isqrts]->fill(n, pow(sph, n));
140      }
141
142      // C & D params
143      const ParisiTensor& parisi = apply<ParisiTensor>(event, "Parisi");
144      const double cparam = parisi.C();
145      const double dparam = parisi.D();
146      _histCParam[_isqrts]->fill(cparam);
147      _histDParam[_isqrts]->fill(dparam);
148      for (int n = 1; n <= 5; ++n) {
149        _histCParamMom[_isqrts]->fill(n, pow(cparam, n));
150      }
151
152      // Hemispheres
153      const Hemispheres& hemi = apply<Hemispheres>(event, "Hemispheres");
154      // The paper says that M_H/L are scaled by sqrt(s), but scaling by E_vis is the way that fits the data...
155      const double hemi_mh = hemi.scaledMhigh();
156      const double hemi_ml = hemi.scaledMlow();
157      /// @todo This shouldn't be necessary... what's going on? Memory corruption suspected :(
158      // if (std::isnan(hemi_ml)) {
159      //   MSG_ERROR("NaN in HemiL! Event = " << numEvents());
160      //   MSG_ERROR(hemi.M2low() << ", " << hemi.E2vis());
161      // }
162      if (!std::isnan(hemi_mh) && !std::isnan(hemi_ml)) {
163        const double hemi_bmax = hemi.Bmax();
164        const double hemi_bmin = hemi.Bmin();
165        const double hemi_bsum = hemi.Bsum();
166        _histHemiMassH[_isqrts]->fill(hemi_mh);
167        _histHemiMassL[_isqrts]->fill(hemi_ml);
168        _histHemiBroadW[_isqrts]->fill(hemi_bmax);
169        _histHemiBroadN[_isqrts]->fill(hemi_bmin);
170        _histHemiBroadT[_isqrts]->fill(hemi_bsum);
171        for (int n = 1; n <= 5; ++n) {
172          // if (std::isnan(pow(hemi_ml, n))) MSG_ERROR("NaN in HemiL moment! Event = " << numEvents());
173          _histHemiMassHMom[_isqrts]->fill(n, pow(hemi_mh, n));
174          _histHemiMassLMom[_isqrts]->fill(n, pow(hemi_ml, n));
175          _histHemiBroadWMom[_isqrts]->fill(n, pow(hemi_bmax, n));
176          _histHemiBroadNMom[_isqrts]->fill(n, pow(hemi_bmin, n));
177          _histHemiBroadTMom[_isqrts]->fill(n, pow(hemi_bsum, n));
178        }
179      }
180    }
181
182
183    void finalize() {
184      scale(_hist1MinusT[_isqrts], 1.0 / *_sumWTrack2);
185      scale(_histTMajor[_isqrts], 1.0 / *_sumWTrack2);
186      scale(_histTMinor[_isqrts], 1.0 / *_sumWTrack2);
187      scale(_histOblateness[_isqrts], 1.0 / *_sumWTrack2);
188      scale(_histSphericity[_isqrts], 1.0 / *_sumWTrack2);
189      scale(_histAplanarity[_isqrts], 1.0 / *_sumWTrack2);
190      scale(_histHemiMassH[_isqrts], 1.0 / *_sumWTrack2);
191      scale(_histHemiMassL[_isqrts], 1.0 / *_sumWTrack2);
192      scale(_histHemiBroadW[_isqrts], 1.0 / *_sumWTrack2);
193      scale(_histHemiBroadN[_isqrts], 1.0 / *_sumWTrack2);
194      scale(_histHemiBroadT[_isqrts], 1.0 / *_sumWTrack2);
195      scale(_histCParam[_isqrts], 1.0 / *_sumWTrack2);
196      scale(_histDParam[_isqrts], 1.0 / *_sumWTrack2);
197      scale(_histY23Durham[_isqrts], 1.0 / *_sumWJet3);
198      //
199      scale(_hist1MinusTMom[_isqrts], 1.0 / *_sumWTrack2);
200      scale(_histTMajorMom[_isqrts], 1.0 / *_sumWTrack2);
201      scale(_histTMinorMom[_isqrts], 1.0 / *_sumWTrack2);
202      scale(_histOblatenessMom[_isqrts], 1.0 / *_sumWTrack2);
203      scale(_histSphericityMom[_isqrts], 1.0 / *_sumWTrack2);
204      scale(_histHemiMassHMom[_isqrts], 1.0 / *_sumWTrack2);
205      scale(_histHemiMassLMom[_isqrts], 1.0 / *_sumWTrack2);
206      scale(_histHemiBroadWMom[_isqrts], 1.0 / *_sumWTrack2);
207      scale(_histHemiBroadNMom[_isqrts], 1.0 / *_sumWTrack2);
208      scale(_histHemiBroadTMom[_isqrts], 1.0 / *_sumWTrack2);
209      scale(_histCParamMom[_isqrts], 1.0 / *_sumWTrack2);
210      scale(_histY23DurhamMom[_isqrts], 1.0 / *_sumWJet3);
211    }
212
213    /// @}
214
215
216  private:
217
218    /// Beam energy index for histograms
219    int _isqrts = -1;
220
221    /// Counters of event weights passing the cuts
222    /// @{
223    CounterPtr _sumWTrack2, _sumWJet3;
224    /// @}
225
226    /// @name Event shape histos at 4 energies
227    /// @{
228    Histo1DPtr _hist1MinusT[4];
229    Histo1DPtr _histHemiMassH[4];
230    Histo1DPtr _histCParam[4];
231    Histo1DPtr _histHemiBroadT[4];
232    Histo1DPtr _histHemiBroadW[4];
233    Histo1DPtr _histY23Durham[4];
234    Histo1DPtr _histTMajor[4];
235    Histo1DPtr _histTMinor[4];
236    Histo1DPtr _histAplanarity[4];
237    Histo1DPtr _histSphericity[4];
238    Histo1DPtr _histOblateness[4];
239    Histo1DPtr _histHemiMassL[4];
240    Histo1DPtr _histHemiBroadN[4];
241    Histo1DPtr _histDParam[4];
242    /// @}
243
244    /// @name Event shape moment histos at 4 energies
245    /// @{
246    Histo1DPtr _hist1MinusTMom[4];
247    Histo1DPtr _histHemiMassHMom[4];
248    Histo1DPtr _histCParamMom[4];
249    Histo1DPtr _histHemiBroadTMom[4];
250    Histo1DPtr _histHemiBroadWMom[4];
251    Histo1DPtr _histY23DurhamMom[4];
252    Histo1DPtr _histTMajorMom[4];
253    Histo1DPtr _histTMinorMom[4];
254    Histo1DPtr _histSphericityMom[4];
255    Histo1DPtr _histOblatenessMom[4];
256    Histo1DPtr _histHemiMassLMom[4];
257    Histo1DPtr _histHemiBroadNMom[4];
258    /// @}
259
260  };
261
262
263
264  RIVET_DECLARE_ALIASED_PLUGIN(OPAL_2004_S6132243, OPAL_2004_I669402);
265
266}