rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

L3_2004_I652683

Jet rates and event shapes at LEP I+II
Experiment: L3 (LEP I+II)
Inspire ID: 652683
Status: VALIDATED
Authors:
  • Adil Jueid
References:
  • Phys.Rept. 399 (2004) 71-174
  • 10.1016/j.physrep.2004.07.002
  • arXiv: hep-ex/0406049
Beams: e+ e-
Beam energies: (20.7, 20.7); (27.6, 27.6); (32.7, 32.7); (37.9, 37.9); (41.1, 41.1); (42.5, 42.5); (45.6, 45.6); (65.0, 65.0); (68.0, 68.0); (80.7, 80.7); (86.2, 86.2); (91.4, 91.4); (94.3, 94.3); (97.2, 97.2); (100.1, 100.1); (103.1, 103.1) GeV
Run details:
  • $e^+ e^- \to q \bar{q}$ at several energies. Beam energy must be specified as analysis option "ENERGY" when rivet-merge'ing samples.

Event shapes, charged particle multiplicity and scaled momentum distributions, flavour separated at the Z-boson peak. The analysis also includes non-flavour separated event shapes, charged particle multiplicity and scaled momentum distributions at events both above and below the Z pole. There are also jet distributions in the paper which are not currently implemented. Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples.

Source code: L3_2004_I652683.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/FinalState.hh"
  6#include "Rivet/Projections/ChargedFinalState.hh"
  7#include "Rivet/Projections/Thrust.hh"
  8#include "Rivet/Projections/ParisiTensor.hh"
  9#include "Rivet/Projections/Hemispheres.hh"
 10
 11#define I_KNOW_THE_INITIAL_QUARKS_PROJECTION_IS_DODGY_BUT_NEED_TO_USE_IT
 12#include "Rivet/Projections/InitialQuarks.hh"
 13#include "fastjet/EECambridgePlugin.hh"
 14
 15namespace Rivet {
 16
 17
 18  /// Jet rates and event shapes at LEP I+II
 19  class L3_2004_I652683 : public Analysis {
 20  public:
 21
 22    /// Constructor
 23    RIVET_DEFAULT_ANALYSIS_CTOR(L3_2004_I652683);
 24
 25    /// Book histograms and initialise projections before the run
 26    void init() {
 27
 28      // Projections to use
 29      const FinalState FS;
 30      declare(FS, "FS");
 31      declare(Beam(), "beams");
 32      const ChargedFinalState CFS;
 33      declare(CFS, "CFS");
 34      const Thrust thrust(FS);
 35      declare(thrust, "thrust");
 36      declare(ParisiTensor(FS), "Parisi");
 37      declare(Hemispheres(thrust), "Hemispheres");
 38      declare(InitialQuarks(), "initialquarks");
 39      FastJets jadeJets = FastJets(FS, FastJets::JADE, 0.7, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
 40      FastJets durhamJets = FastJets(FS, FastJets::DURHAM, 0.7, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
 41      FastJets cambridgeJets = FastJets(FS, FastJets::CAM, 0.7, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
 42      declare(jadeJets, "JadeJets");
 43      declare(durhamJets, "DurhamJets");
 44
 45      // Book the histograms
 46      if (isCompatibleWithSqrtS(91.2)) {
 47        // z pole
 48        book(_h_Thrust_udsc             , 47, 1, 1);
 49        book(_h_Thrust_bottom           , 47, 1, 2);
 50        book(_h_heavyJetmass_udsc       , 48, 1, 1);
 51        book(_h_heavyJetmass_bottom     , 48, 1, 2);
 52        book(_h_totalJetbroad_udsc      , 49, 1, 1);
 53        book(_h_totalJetbroad_bottom    , 49, 1, 2);
 54        book(_h_wideJetbroad_udsc       , 50, 1, 1);
 55        book(_h_wideJetbroad_bottom     , 50, 1, 2);
 56        book(_h_Cparameter_udsc         , 51, 1, 1);
 57        book(_h_Cparameter_bottom       , 51, 1, 2);
 58        book(_h_Dparameter_udsc         , 52, 1, 1);
 59        book(_h_Dparameter_bottom       , 52, 1, 2);
 60        book(_h_Ncharged                , "/TMP/NCHARGED"     , 28, 1, 57);
 61        book(_h_Ncharged_udsc           , "/TMP/NCHARGED_UDSC", 28, 1, 57);
 62        book(_h_Ncharged_bottom         , "/TMP/NCHARGED_B"   , 27, 3, 57);
 63        book(_h_scaledMomentum          , 65, 1, 1);
 64        book(_h_scaledMomentum_udsc     , 65, 1, 2);
 65        book(_h_scaledMomentum_bottom   ,  65, 1, 3);
 66      }
 67      else if(sqrtS()/GeV < 90) {
 68        int i1(-1),i2(-1);
 69        if(isCompatibleWithSqrtS(41.4)) {
 70          i1=0; i2=1;
 71        }
 72        else if(isCompatibleWithSqrtS(55.3)) {
 73          i1=0; i2=2;
 74        }
 75        else if(isCompatibleWithSqrtS(65.4)) {
 76          i1=0; i2=3;
 77        }
 78        else if(isCompatibleWithSqrtS(75.7)) {
 79          i1=1; i2=1;
 80        }
 81        else if(isCompatibleWithSqrtS(82.3)) {
 82          i1=1; i2=2;
 83        }
 84        else if(isCompatibleWithSqrtS(85.1)) {
 85          i1=1; i2=3;
 86        }
 87        else
 88          MSG_ERROR("Beam energy not supported!");
 89        book(_h_thrust , 21+i1,1,i2);
 90        book(_h_rho    , 26+i1,1,i2);
 91        book(_h_B_T    , 31+i1,1,i2);
 92        book(_h_B_W    , 36+i1,1,i2);
 93      }
 94      else if(sqrtS()/GeV > 120) {
 95        int i1(-1),i2(-1);
 96        if(isCompatibleWithSqrtS(130.1)) {
 97          i1=0; i2=1;
 98        }
 99        else if(isCompatibleWithSqrtS(136.1)) {
100          i1=0; i2=2;
101        }
102        else if(isCompatibleWithSqrtS(161.3)) {
103          i1=0; i2=3;
104        }
105        else if(isCompatibleWithSqrtS(172.3)) {
106          i1=1; i2=1;
107        }
108        else if(isCompatibleWithSqrtS(182.8)) {
109          i1=1; i2=2;
110        }
111        else if(isCompatibleWithSqrtS(188.6)) {
112          i1=1; i2=3;
113        }
114        else if(isCompatibleWithSqrtS(194.4)) {
115          i1=2; i2=1;
116        }
117        else if(isCompatibleWithSqrtS(200.2)) {
118          i1=2; i2=2;
119        }
120        else if(isCompatibleWithSqrtS(206.2)) {
121          i1=2; i2=3;
122        }
123        else
124          MSG_ERROR("Beam energy not supported!");
125	book(_h_thrust , 23+i1,1,i2);
126	book(_h_rho    , 28+i1,1,i2);
127	book(_h_B_T    , 33+i1,1,i2);
128	book(_h_B_W    , 38+i1,1,i2);
129	book(_h_C      , 41+i1,1,i2);
130	book(_h_D      , 44+i1,1,i2);
131	book(_h_N      , "/TMP/NCHARGED", 22, 9, 53);
132	book(_h_xi     , 66+i1,1,i2);
133	// and the jets
134	int i3 = 3*i1+i2;
135	book(_h_y_2_JADE  ,   i3,1,1);
136	book(_h_y_3_JADE  ,   i3,1,2);
137	book(_h_y_4_JADE  ,   i3,1,3);
138	book(_h_y_5_JADE  ,   i3,1,4);
139	book(_h_y_2_Durham, 9+i3,1,1);
140	book(_h_y_3_Durham, 9+i3,1,2);
141	book(_h_y_4_Durham, 9+i3,1,3);
142	book(_h_y_5_Durham, 9+i3,1,4);
143	if(i3==8||i3==9) {
144	  book(_h_y_2_Cambridge,11+i3,1,1);
145	  book(_h_y_3_Cambridge,11+i3,1,2);
146	  book(_h_y_4_Cambridge,11+i3,1,3);
147	  book(_h_y_5_Cambridge,11+i3,1,4);
148	}
149      }
150
151      book(_sumW_udsc, "_sumW_udsc");
152      book(_sumW_b, "_sumW_b");
153      book(_sumW_ch, "_sumW_ch");
154      book(_sumW_ch_udsc, "_sumW_ch_udsc");
155      book(_sumW_ch_b, "_sumW_ch_b");
156
157    }
158
159
160    /// Perform the per-event analysis
161    void analyze(const Event& event) {
162      // Get beam average momentum
163      const ParticlePair& beams = apply<Beam>(event, "beams").beams();
164      const double beamMomentum = ( beams.first.p3().mod() + beams.second.p3().mod() ) / 2.0;
165
166      // InitialQuarks projection to have udsc events separated from b events
167      /// @todo Yuck!!! Eliminate when possible...
168      int iflav = 0;
169      // only need the flavour at Z pole
170      if(_h_Thrust_udsc) {
171	int flavour = 0;
172	const InitialQuarks& iqf = apply<InitialQuarks>(event, "initialquarks");
173	Particles quarks;
174	if ( iqf.particles().size() == 2 ) {
175	  flavour = iqf.particles().front().abspid();
176	  quarks  = iqf.particles();
177	} else {
178	  map<int, Particle> quarkmap;
179	  for (const Particle& p : iqf.particles()) {
180	    if (quarkmap.find(p.pid()) == quarkmap.end()) quarkmap[p.pid()] = p;
181	    else if (quarkmap[p.pid()].E() < p.E()) quarkmap[p.pid()] = p;
182	  }
183	  double max_energy = 0.;
184	  for (int i = 1; i <= 5; ++i) {
185	    double energy = 0.;
186	    if (quarkmap.find(i) != quarkmap.end())
187	      energy += quarkmap[ i].E();
188	    if (quarkmap.find(-i) != quarkmap.end())
189	      energy += quarkmap[-i].E();
190	    if (energy > max_energy)
191	      flavour = i;
192	  }
193	  if (quarkmap.find(flavour) != quarkmap.end())
194	    quarks.push_back(quarkmap[flavour]);
195	  if (quarkmap.find(-flavour) != quarkmap.end())
196	    quarks.push_back(quarkmap[-flavour]);
197	}
198	// Flavour label
199	/// @todo Change to a bool?
200	iflav = (flavour == PID::DQUARK || flavour == PID::UQUARK || flavour == PID::SQUARK || flavour == PID::CQUARK) ? 1 : (flavour == PID::BQUARK) ? 5 : 0;
201      }
202      // Update weight sums
203      if (iflav == 1) {
204        _sumW_udsc->fill();
205      } else if (iflav == 5) {
206        _sumW_b->fill();
207      }
208      _sumW_ch->fill();
209
210      // Charged multiplicity
211      const FinalState& cfs = applyProjection<FinalState>(event, "CFS");
212      if(_h_Ncharged) _h_Ncharged->fill(cfs.size());
213      if (iflav == 1) {
214        _sumW_ch_udsc->fill();
215        _h_Ncharged_udsc->fill(cfs.size());
216      } else if (iflav == 5) {
217        _sumW_ch_b->fill();
218        _h_Ncharged_bottom->fill(cfs.size());
219      }
220      else if(_h_N) {
221	_h_N->fill(cfs.size());
222      }
223
224      // Scaled momentum
225      const Particles& chparticles = cfs.particlesByPt();
226      for (const Particle& p : chparticles) {
227        const Vector3 momentum3 = p.p3();
228        const double mom = momentum3.mod();
229        const double scaledMom = mom/beamMomentum;
230        const double logScaledMom = std::log(scaledMom);
231        if(_h_scaledMomentum) _h_scaledMomentum->fill(-logScaledMom);
232        if (iflav == 1) {
233          _h_scaledMomentum_udsc->fill(-logScaledMom);
234        } else if (iflav == 5) {
235          _h_scaledMomentum_bottom->fill(-logScaledMom);
236        }
237	else if(_h_xi) {
238	  _h_xi->fill(-logScaledMom);
239	}
240      }
241
242      // Thrust
243      const Thrust& thrust = applyProjection<Thrust>(event, "thrust");
244      if (iflav == 1) {
245        _h_Thrust_udsc->fill(thrust.thrust());
246      } else if (iflav == 5) {
247        _h_Thrust_bottom->fill(thrust.thrust());
248      }
249      else if(_h_thrust) {
250        _h_thrust->fill(1.-thrust.thrust());
251      }
252
253      // C and D Parisi parameters
254      const ParisiTensor& parisi = applyProjection<ParisiTensor>(event, "Parisi");
255      if (iflav == 1) {
256        _h_Cparameter_udsc->fill(parisi.C());
257        _h_Dparameter_udsc->fill(parisi.D());
258      } else if (iflav == 5) {
259        _h_Cparameter_bottom->fill(parisi.C());
260        _h_Dparameter_bottom->fill(parisi.D());
261      }
262      else if(_h_C) {
263	_h_C->fill(parisi.C());
264	_h_D->fill(parisi.D());
265      }
266
267      // The hemisphere variables
268      const Hemispheres& hemisphere = applyProjection<Hemispheres>(event, "Hemispheres");
269      if (iflav == 1) {
270        _h_heavyJetmass_udsc->fill(hemisphere.scaledM2high());
271        _h_totalJetbroad_udsc->fill(hemisphere.Bsum());
272        _h_wideJetbroad_udsc->fill(hemisphere.Bmax());
273      } else if (iflav == 5) {
274        _h_heavyJetmass_bottom->fill(hemisphere.scaledM2high());
275        _h_totalJetbroad_bottom->fill(hemisphere.Bsum());
276        _h_wideJetbroad_bottom->fill(hemisphere.Bmax());
277      }
278      else if (_h_rho) {
279        _h_rho->fill(hemisphere.scaledM2high());
280        _h_B_T->fill(hemisphere.Bsum());
281        _h_B_W->fill(hemisphere.Bmax());
282      }
283      // jade jet rates
284      if(_h_y_2_JADE) {
285	const FastJets& jadejet = apply<FastJets>(event, "JadeJets");
286	if (jadejet.clusterSeq()) {
287	  const double y_23 = jadejet.clusterSeq()->exclusive_ymerge_max(2);
288	  const double y_34 = jadejet.clusterSeq()->exclusive_ymerge_max(3);
289	  const double y_45 = jadejet.clusterSeq()->exclusive_ymerge_max(4);
290	  const double y_56 = jadejet.clusterSeq()->exclusive_ymerge_max(5);
291	  for (size_t i = 0; i < _h_y_2_JADE->numBins(); ++i) {
292	    double ycut = _h_y_2_JADE->bin(i).xMid();
293	    double width = _h_y_2_JADE->bin(i).width();
294	    if (y_23 < ycut) _h_y_2_JADE->fillBin(i,width);
295	  }
296	  for (size_t i = 0; i < _h_y_3_JADE->numBins(); ++i) {
297	    double ycut = _h_y_3_JADE->bin(i).xMid();
298	    double width = _h_y_3_JADE->bin(i).width();
299	    if (y_34 < ycut && y_23 > ycut) {
300	      _h_y_3_JADE->fillBin(i,width);
301	    }
302	  }
303	  for (size_t i = 0; i < _h_y_4_JADE->numBins(); ++i) {
304	    double ycut = _h_y_4_JADE->bin(i).xMid();
305	    double width = _h_y_4_JADE->bin(i).width();
306	    if (y_45 < ycut && y_34 > ycut) {
307	      _h_y_4_JADE->fillBin(i,width);
308	    }
309	  }
310	  for (size_t i = 0; i < _h_y_5_JADE->numBins(); ++i) {
311	    double ycut = _h_y_5_JADE->bin(i).xMid();
312	    double width = _h_y_5_JADE->bin(i).width();
313	    if (y_56 < ycut && y_45 > ycut) {
314	      _h_y_5_JADE->fillBin(i,width);
315	    }
316	  }
317	}
318      }
319      // Durham jet rates
320      if(_h_y_2_Durham) {
321	const FastJets& durhamjet = apply<FastJets>(event, "DurhamJets");
322	if (durhamjet.clusterSeq()) {
323	  const double y_23 = durhamjet.clusterSeq()->exclusive_ymerge_max(2);
324	  const double y_34 = durhamjet.clusterSeq()->exclusive_ymerge_max(3);
325	  const double y_45 = durhamjet.clusterSeq()->exclusive_ymerge_max(4);
326	  const double y_56 = durhamjet.clusterSeq()->exclusive_ymerge_max(5);
327	  for (size_t i = 0; i < _h_y_2_Durham->numBins(); ++i) {
328	    double ycut = _h_y_2_Durham->bin(i).xMid();
329	    double width = _h_y_2_Durham->bin(i).width();
330	    if (y_23 < ycut) _h_y_2_Durham->fillBin(i,width);
331	  }
332	  for (size_t i = 0; i < _h_y_3_Durham->numBins(); ++i) {
333	    double ycut = _h_y_3_Durham->bin(i).xMid();
334	    double width = _h_y_3_Durham->bin(i).width();
335	    if (y_34 < ycut && y_23 > ycut) {
336	      _h_y_3_Durham->fillBin(i,width);
337	    }
338	  }
339	  for (size_t i = 0; i < _h_y_4_Durham->numBins(); ++i) {
340	    double ycut = _h_y_4_Durham->bin(i).xMid();
341	    double width = _h_y_4_Durham->bin(i).width();
342	    if (y_45 < ycut && y_34 > ycut) {
343	      _h_y_4_Durham->fillBin(i,width);
344	    }
345	  }
346	  for (size_t i = 0; i < _h_y_5_Durham->numBins(); ++i) {
347	    double ycut = _h_y_5_Durham->bin(i).xMid();
348	    double width = _h_y_5_Durham->bin(i).width();
349	    if (y_56 < ycut && y_45 > ycut) {
350	      _h_y_5_Durham->fillBin(i,width);
351	    }
352	  }
353	}
354      }
355      // Cambridge
356      if(_h_y_2_Cambridge) {
357	PseudoJets pjs;
358	const FinalState& fs = applyProjection<FinalState>(event, "FS");
359	for (size_t i = 0; i < fs.particles().size(); ++i) {
360	  fastjet::PseudoJet pj = fs.particles()[i];
361	  pjs.push_back(pj);
362	}
363	for (size_t i = 0; i < _h_y_2_Cambridge->numBins(); ++i) {
364	  double ycut = _h_y_2_Cambridge->bin(i).xMid();
365	  double width = _h_y_2_Cambridge->bin(i).width();
366	  fastjet::EECambridgePlugin plugin(ycut);
367	  fastjet::JetDefinition jdef(&plugin);
368	  fastjet::ClusterSequence cseq(pjs, jdef);
369	  unsigned int njet = cseq.inclusive_jets().size();
370	  if(njet==2)
371	    _h_y_2_Cambridge->fillBin(i,width);
372	  else if(njet==3) {
373	    if(i<_h_y_3_Cambridge->numBins()) _h_y_3_Cambridge->fillBin(i,width);
374	  }
375	  else if(njet==4) {
376	    if(i<_h_y_4_Cambridge->numBins()) _h_y_4_Cambridge->fillBin(i,width);
377	  }
378	  else if(njet==5) {
379	    if(i<_h_y_5_Cambridge->numBins()) _h_y_5_Cambridge->fillBin(i,width);
380	  }
381	}
382      }
383    }
384    
385    Scatter2DPtr convertHisto(unsigned int ix,unsigned int iy, unsigned int iz, Histo1DPtr histo) {
386      Scatter2D temphisto(refData(ix, iy, iz));
387      Scatter2DPtr mult;
388      book(mult, ix, iy, iz);
389      for (size_t b = 0; b < temphisto.numPoints(); b++) {
390	const double x  = temphisto.point(b).x();
391	pair<double,double> ex = temphisto.point(b).xErrs();
392	double y    = histo->bins()[b].area();
393	double yerr = histo->bins()[b].areaErr();
394	mult->addPoint(x, y, ex, make_pair(yerr,yerr));
395      }
396      return mult;
397    }
398
399    /// Normalise histograms etc., after the run
400    void finalize() {
401      // Z pole plots
402      if(_h_Thrust_udsc) {
403        scale(_h_Thrust_udsc,  1/_sumW_udsc->sumW());
404        scale(_h_heavyJetmass_udsc,  1/_sumW_udsc->sumW());
405        scale(_h_totalJetbroad_udsc, 1/_sumW_udsc->sumW());
406	scale(_h_wideJetbroad_udsc, 1/_sumW_udsc->sumW());
407        scale(_h_Cparameter_udsc, 1/_sumW_udsc->sumW());
408        scale(_h_Dparameter_udsc, 1/_sumW_udsc->sumW());
409        scale(_h_Thrust_bottom, 1./_sumW_b->sumW());
410        scale(_h_heavyJetmass_bottom, 1./_sumW_b->sumW());
411        scale(_h_totalJetbroad_bottom, 1./_sumW_b->sumW());
412	scale(_h_wideJetbroad_bottom, 1./_sumW_b->sumW());
413        scale(_h_Cparameter_bottom, 1./_sumW_b->sumW());
414        scale(_h_Dparameter_bottom, 1./_sumW_b->sumW());
415        scale(_h_Ncharged, 1./_sumW_ch->sumW());
416        scale(_h_Ncharged_udsc, 1./_sumW_ch_udsc->sumW());
417        scale(_h_Ncharged_bottom, 1./_sumW_ch_b->sumW());
418        convertHisto(59,1,1,_h_Ncharged       );
419        convertHisto(59,1,2,_h_Ncharged_udsc  );
420        convertHisto(59,1,3,_h_Ncharged_bottom);
421
422        scale(_h_scaledMomentum, 1/_sumW_ch->sumW());
423        scale(_h_scaledMomentum_udsc, 1/_sumW_ch_udsc->sumW());
424        scale(_h_scaledMomentum_bottom, 1/_sumW_ch_b->sumW());
425      }
426      else {
427	if(_h_thrust) normalize(_h_thrust);
428	if(_h_rho) normalize(_h_rho);
429	if(_h_B_T) normalize(_h_B_T);
430	if(_h_B_W) normalize(_h_B_W);
431	if(_h_C) normalize(_h_C);
432	if(_h_D) normalize(_h_D);
433	if(_h_N) normalize(_h_N);
434	if(_h_xi) scale(_h_xi,1./sumOfWeights());
435	
436      
437      Scatter2DPtr mult;
438	if(_h_N) {
439    if(isCompatibleWithSqrtS(130.1)) {
440	    convertHisto(60, 1, 1, _h_N);
441	  }
442	  else if(isCompatibleWithSqrtS(136.1)) {
443	    convertHisto(60, 1, 2, _h_N);
444	  }
445	  else if(isCompatibleWithSqrtS(161.3)) {
446	    convertHisto(60, 1, 3, _h_N);
447	  }
448	  else if(isCompatibleWithSqrtS(172.3)) {
449	    convertHisto(61, 1, 1, _h_N);
450	  }
451	  else if(isCompatibleWithSqrtS(182.8)) {
452	    convertHisto(61, 1, 2, _h_N);
453	  }
454	  else if(isCompatibleWithSqrtS(188.6)) {
455	    convertHisto(61, 1, 3, _h_N);
456	  }
457	  else if(isCompatibleWithSqrtS(194.4)) {
458	    convertHisto(62, 1, 1, _h_N);
459	  }
460	  else if(isCompatibleWithSqrtS(200.2)) {
461	    convertHisto(62, 1, 2, _h_N);
462	  }
463	  else if(isCompatibleWithSqrtS(206.2)) {
464	    convertHisto(62, 1, 3, _h_N);
465	  }
466	}
467	// // the jets
468	if(_h_y_2_JADE) {
469	  scale(_h_y_2_JADE, 1./ sumOfWeights());
470	  scale(_h_y_3_JADE, 1./ sumOfWeights());
471	  scale(_h_y_4_JADE, 1./ sumOfWeights());
472	  scale(_h_y_5_JADE, 1./ sumOfWeights());
473	}
474	if(_h_y_2_Durham) {
475	  scale(_h_y_2_Durham, 1./ sumOfWeights());
476	  scale(_h_y_3_Durham, 1./ sumOfWeights());
477	  scale(_h_y_4_Durham, 1./ sumOfWeights());
478	  scale(_h_y_5_Durham, 1./ sumOfWeights());
479	}
480	if(_h_y_2_Cambridge) {
481	  scale(_h_y_2_Cambridge, 1./ sumOfWeights());
482	  scale(_h_y_3_Cambridge, 1./ sumOfWeights());
483	  scale(_h_y_4_Cambridge, 1./ sumOfWeights());
484	  scale(_h_y_5_Cambridge, 1./ sumOfWeights());
485	}
486      }
487    }
488
489
490    /// Weight counters
491    CounterPtr _sumW_udsc, _sumW_b, _sumW_ch, _sumW_ch_udsc, _sumW_ch_b;
492
493    /// @name Histograms
494    //@{
495    // at the Z pole
496    Histo1DPtr _h_Thrust_udsc, _h_Thrust_bottom;
497    Histo1DPtr _h_heavyJetmass_udsc, _h_heavyJetmass_bottom;
498    Histo1DPtr _h_totalJetbroad_udsc, _h_totalJetbroad_bottom;
499    Histo1DPtr _h_wideJetbroad_udsc, _h_wideJetbroad_bottom;
500    Histo1DPtr _h_Cparameter_udsc, _h_Cparameter_bottom;
501    Histo1DPtr _h_Dparameter_udsc, _h_Dparameter_bottom;
502    Histo1DPtr _h_Ncharged, _h_Ncharged_udsc, _h_Ncharged_bottom;
503    Histo1DPtr _h_scaledMomentum, _h_scaledMomentum_udsc, _h_scaledMomentum_bottom;
504    // at other enegies
505    Histo1DPtr _h_thrust,_h_rho,_h_B_T,_h_B_W,_h_C,_h_D,_h_N,_h_xi;
506    // and the jets
507    Histo1DPtr _h_y_2_JADE,_h_y_3_JADE,_h_y_4_JADE,_h_y_5_JADE;
508    Histo1DPtr _h_y_2_Durham,_h_y_3_Durham,_h_y_4_Durham,_h_y_5_Durham;
509    Histo1DPtr _h_y_2_Cambridge,_h_y_3_Cambridge,_h_y_4_Cambridge,_h_y_5_Cambridge;
510    //@}
511
512  };
513
514
515  // The hook for the plugin system
516  RIVET_DECLARE_PLUGIN(L3_2004_I652683);
517
518}