FastJet  3.3.0
JetDefinition.cc
1 //FJSTARTHEADER
2 // $Id: JetDefinition.cc 4074 2016-03-08 09:09:25Z soyez $
3 //
4 // Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5 //
6 //----------------------------------------------------------------------
7 // This file is part of FastJet.
8 //
9 // FastJet is free software; you can redistribute it and/or modify
10 // it under the terms of the GNU General Public License as published by
11 // the Free Software Foundation; either version 2 of the License, or
12 // (at your option) any later version.
13 //
14 // The algorithms that underlie FastJet have required considerable
15 // development. They are described in the original FastJet paper,
16 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17 // FastJet as part of work towards a scientific publication, please
18 // quote the version you use and include a citation to the manual and
19 // optionally also to hep-ph/0512210.
20 //
21 // FastJet is distributed in the hope that it will be useful,
22 // but WITHOUT ANY WARRANTY; without even the implied warranty of
23 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 // GNU General Public License for more details.
25 //
26 // You should have received a copy of the GNU General Public License
27 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28 //----------------------------------------------------------------------
29 //FJENDHEADER
30 
31 #include "fastjet/JetDefinition.hh"
32 #include "fastjet/Error.hh"
33 #include "fastjet/CompositeJetStructure.hh"
34 #include<sstream>
35 
36 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
37 
38 using namespace std;
39 
40 const double JetDefinition::max_allowable_R = 1000.0;
41 
42 //----------------------------------------------------------------------
43 // [NB: implementation was getting complex, so in 2.4-devel moved it
44 // from .hh to .cc]
45 JetDefinition::JetDefinition(JetAlgorithm jet_algorithm_in,
46  double R_in,
47  RecombinationScheme recomb_scheme_in,
48  Strategy strategy_in,
49  int nparameters) :
50  _jet_algorithm(jet_algorithm_in), _Rparam(R_in), _strategy(strategy_in) {
51 
52  // set R parameter or ensure its sensibleness, as appropriate
53  if (_jet_algorithm == ee_kt_algorithm) {
54  _Rparam = 4.0; // introduce a fictional R that ensures that
55  // our clustering sequence will not produce
56  // "beam" jets except when only a single particle remains.
57  // Any value > 2 would have done here
58  } else {
59  // We maintain some limit on R because particles with pt=0, m=0
60  // can have rapidities O(100000) and one doesn't want the
61  // clustering to start including them as if their rapidities were
62  // physical.
63  if (R_in > max_allowable_R) {
64  ostringstream oss;
65  oss << "Requested R = " << R_in << " for jet definition is larger than max_allowable_R = " << max_allowable_R;
66  throw Error(oss.str());
67  }
68  }
69 
70  // cross-check the number of parameters that were declared in setting up the
71  // algorithm (passed internally from the public constructors)
72  unsigned int nparameters_expected = n_parameters_for_algorithm(jet_algorithm_in);
73  if (nparameters != (int) nparameters_expected){
74  ostringstream oss;
75  oss << "The jet algorithm you requested ("
76  << jet_algorithm_in << ") should be constructed with " << nparameters_expected
77  << " parameter(s) but was called with " << nparameters << " parameter(s)\n";
78  throw Error(oss.str());
79  }
80 
81  // make sure the strategy requested is sensible
82  assert (_strategy != plugin_strategy);
83 
84  _plugin = NULL;
85  set_recombination_scheme(recomb_scheme_in);
86  set_extra_param(0.0); // make sure it's defined
87 }
88 
89 
90 //----------------------------------------------------------------------
91 // returns true if the jet definition involves an algorithm
92 // intended for use on a spherical geometry (e.g. e+e- algorithms,
93 // as opposed to most pp algorithms, which use a cylindrical,
94 // rapidity-phi geometry).
97  return plugin()->is_spherical();
98  } else {
99  return (jet_algorithm() == ee_kt_algorithm || // as of 2013-02-14, the two
100  jet_algorithm() == ee_genkt_algorithm // native spherical algorithms
101  );
102  }
103 }
104 
105 //----------------------------------------------------------------------
107  ostringstream name;
108 
109  name << description_no_recombiner();
110 
112  return name.str();
113  }
114 
116  name << " with ";
117  else
118  name << " and ";
119  name << recombiner()->description();
120 
121  return name.str();
122 }
123 
124 //----------------------------------------------------------------------
126 
127  ostringstream name;
128  if (jet_algorithm() == plugin_algorithm) {
129  return plugin()->description();
130  } else if (jet_algorithm() == undefined_jet_algorithm) {
131  return "uninitialised JetDefinition (jet_algorithm=undefined_jet_algorithm)" ;
132  }
133 
136  case 0: name << " (NB: no R)"; break;
137  case 1: name << " with R = " << R(); break; // the parameter is always R
138  case 2:
139  // the 1st parameter is always R
140  name << " with R = " << R();
141  // the 2nd depends on the algorithm
143  name << "and a special hack whereby particles with kt < "
144  << extra_param() << "are treated as passive ghosts";
145  } else {
146  name << ", p = " << extra_param();
147  }
148  };
149 
150  return name.str();
151 }
152 
153 //----------------------------------------------------------------------
155  ostringstream name;
156  switch (jet_alg){
157  case plugin_algorithm: return "plugin algorithm";
158  case kt_algorithm: return "Longitudinally invariant kt algorithm";
159  case cambridge_algorithm: return "Longitudinally invariant Cambridge/Aachen algorithm";
160  case antikt_algorithm: return "Longitudinally invariant anti-kt algorithm";
161  case genkt_algorithm: return "Longitudinally invariant generalised kt algorithm";
162  case cambridge_for_passive_algorithm: return "Longitudinally invariant Cambridge/Aachen algorithm";
163  case ee_kt_algorithm: return "e+e- kt (Durham) algorithm (NB: no R)";
164  case ee_genkt_algorithm: return "e+e- generalised kt algorithm";
165  case undefined_jet_algorithm: return "undefined jet algorithm";
166  default:
167  throw Error("JetDefinition::algorithm_description(): unrecognized jet_algorithm");
168  };
169 }
170 
171 //----------------------------------------------------------------------
173  switch (jet_alg) {
174  case ee_kt_algorithm: return 0;
175  case genkt_algorithm:
176  case ee_genkt_algorithm: return 2;
177  default: return 1;
178  };
179 }
180 
181 //----------------------------------------------------------------------
183  RecombinationScheme recomb_scheme) {
184  _default_recombiner = JetDefinition::DefaultRecombiner(recomb_scheme);
185 
186  // do not forget to delete the existing recombiner if needed
187  if (_shared_recombiner) _shared_recombiner.reset();
188 
189  _recombiner = 0;
190 }
191 
192 void JetDefinition::set_recombiner(const JetDefinition &other_jet_def){
193  // make sure the "invariants" of the other jet def are sensible
194  assert(other_jet_def._recombiner ||
195  other_jet_def.recombination_scheme() != external_scheme);
196 
197  // first treat the situation where we're using the default recombiner
198  if (other_jet_def._recombiner == 0){
199  set_recombination_scheme(other_jet_def.recombination_scheme());
200  return;
201  }
202 
203  // in other cases, copy the pointer to the recombiner
204  _recombiner = other_jet_def._recombiner;
205  // set the default recombiner appropriately
206  _default_recombiner = DefaultRecombiner(external_scheme);
207  // and set the _shared_recombiner to the same state
208  // as in the other_jet_def, whatever that was
209  _shared_recombiner.reset(other_jet_def._shared_recombiner);
210 
211  // NB: it is tempting to go via set_recombiner and then to sort
212  // out the shared part, but this would be dangerous in the
213  // specific (rare?) case where other_jet_def is the same as this
214  // it deletes_recombiner_when_unused. In that case the shared
215  // pointer reset would delete the recombiner.
216 }
217 
218 
219 // returns true if the current jet definitions shares the same
220 // recombiner as teh one passed as an argument
222  // first make sure that they have the same recombination scheme
223  const RecombinationScheme & scheme = recombination_scheme();
224  if (other_jd.recombination_scheme() != scheme) return false;
225 
226  // if the scheme is "external", also check that they have the same
227  // recombiner
228  return (scheme != external_scheme)
229  || (recombiner() == other_jd.recombiner());
230 }
231 
232 /// causes the JetDefinition to handle the deletion of the
233 /// recombiner when it is no longer used
235  if (_recombiner == 0){
236  throw Error("tried to call JetDefinition::delete_recombiner_when_unused() for a JetDefinition without a user-defined recombination scheme");
237  } else if (_shared_recombiner.get()) {
238  throw Error("Error in JetDefinition::delete_recombiner_when_unused: the recombiner is already scheduled for deletion when unused (or was already set as shared)");
239  }
240 
241  _shared_recombiner.reset(_recombiner);
242 }
243 
244 /// allows to let the JetDefinition handle the deletion of the
245 /// plugin when it is no longer used
247  if (_plugin == 0){
248  throw Error("tried to call JetDefinition::delete_plugin_when_unused() for a JetDefinition without a plugin");
249  }
250 
251  _plugin_shared.reset(_plugin);
252 }
253 
254 
255 
257  switch(_recomb_scheme) {
258  case E_scheme:
259  return "E scheme recombination";
260  case pt_scheme:
261  return "pt scheme recombination";
262  case pt2_scheme:
263  return "pt2 scheme recombination";
264  case Et_scheme:
265  return "Et scheme recombination";
266  case Et2_scheme:
267  return "Et2 scheme recombination";
268  case BIpt_scheme:
269  return "boost-invariant pt scheme recombination";
270  case BIpt2_scheme:
271  return "boost-invariant pt2 scheme recombination";
272  case WTA_pt_scheme:
273  return "pt-ordered Winner-Takes-All recombination";
274  // Energy-ordering can lead to dangerous situations with particles at
275  // rest. We instead implement the WTA_modp_scheme
276  //
277  // case WTA_E_scheme:
278  // return "energy-ordered Winner-Takes-All recombination";
279  case WTA_modp_scheme:
280  return "|3-momentum|-ordered Winner-Takes-All recombination";
281  default:
282  ostringstream err;
283  err << "DefaultRecombiner: unrecognized recombination scheme "
284  << _recomb_scheme;
285  throw Error(err.str());
286  }
287 }
288 
289 
291  const PseudoJet & pa, const PseudoJet & pb,
292  PseudoJet & pab) const {
293 
294  double weighta, weightb;
295 
296  switch(_recomb_scheme) {
297  case E_scheme:
298  // a call to reset turns out to be somewhat more efficient
299  // than a sum and assignment
300  //pab = pa + pb;
301  pab.reset(pa.px()+pb.px(),
302  pa.py()+pb.py(),
303  pa.pz()+pb.pz(),
304  pa.E ()+pb.E ());
305  return;
306  // all remaining schemes are massless recombinations and locally
307  // we just set weights, while the hard work is done below...
308  case pt_scheme:
309  case Et_scheme:
310  case BIpt_scheme:
311  weighta = pa.perp();
312  weightb = pb.perp();
313  break;
314  case pt2_scheme:
315  case Et2_scheme:
316  case BIpt2_scheme:
317  weighta = pa.perp2();
318  weightb = pb.perp2();
319  break;
320  case WTA_pt_scheme:{
321  const PseudoJet & phard = (pa.pt2() >= pb.pt2()) ? pa : pb;
322  /// keep y,phi and m from the hardest, sum pt
323  pab.reset_PtYPhiM(pa.pt()+pb.pt(),
324  phard.rap(), phard.phi(), phard.m());
325  return;}
326  // Energy-ordering can lead to dangerous situations with particles at
327  // rest. We instead implement the WTA_modp_scheme
328  //
329  // case WTA_E_scheme:{
330  // const PseudoJet & phard = (pa.E() >= pb.E()) ? pa : pb;
331  // /// keep 3-momentum direction and mass from the hardest, sum energies
332  // ///
333  // /// If the particle with the largest energy is at rest, the sum
334  // /// remains at rest, implying that the mass of the sum is larger
335  // /// than the mass of pa.
336  // double Eab = pa.E() + pb.E();
337  // double scale = (phard.modp2()==0.0)
338  // ? 0.0
339  // : sqrt((Eab*Eab - phard.m2())/phard.modp2());
340  // pab.reset(phard.px()*scale, phard.py()*scale, phard.pz()*scale, Eab);
341  // return;}
342  case WTA_modp_scheme:{
343  // Note: we need to compute both a and b modp. And we need pthard
344  // and its modp. If we want to avoid repeating the test and do
345  // only 2 modp calculations, we'd have to duplicate the code (or
346  // use a pair<const PJ&>). An alternative is to write modp_soft as
347  // modp_ab-modp_hard but this could suffer from larger rounding
348  // errors
349  bool a_hardest = (pa.modp2() >= pb.modp2());
350  const PseudoJet & phard = a_hardest ? pa : pb;
351  const PseudoJet & psoft = a_hardest ? pb : pa;
352  /// keep 3-momentum direction and mass from the hardest, sum modp
353  ///
354  /// If the hardest particle is at rest, the sum remains at rest
355  /// (the energy of the sum is therefore the mass of pa)
356  double modp_hard = phard.modp();
357  double modp_ab = modp_hard + psoft.modp();
358  if (phard.modp2()==0.0){
359  pab.reset(0.0, 0.0, 0.0, phard.m());
360  } else {
361  double scale = modp_ab/modp_hard;
362  pab.reset(phard.px()*scale, phard.py()*scale, phard.pz()*scale,
363  sqrt(modp_ab*modp_ab + phard.m2()));
364  }
365  return;}
366  default:
367  ostringstream err;
368  err << "DefaultRecombiner: unrecognized recombination scheme "
369  << _recomb_scheme;
370  throw Error(err.str());
371  }
372 
373  double perp_ab = pa.perp() + pb.perp();
374  if (perp_ab != 0.0) { // weights also non-zero...
375  double y_ab = (weighta * pa.rap() + weightb * pb.rap())/(weighta+weightb);
376 
377  // take care with periodicity in phi...
378  double phi_a = pa.phi(), phi_b = pb.phi();
379  if (phi_a - phi_b > pi) phi_b += twopi;
380  if (phi_a - phi_b < -pi) phi_b -= twopi;
381  double phi_ab = (weighta * phi_a + weightb * phi_b)/(weighta+weightb);
382 
383  // this is much more efficient...
384  pab.reset_PtYPhiM(perp_ab,y_ab,phi_ab);
385  // pab = PseudoJet(perp_ab*cos(phi_ab),
386  // perp_ab*sin(phi_ab),
387  // perp_ab*sinh(y_ab),
388  // perp_ab*cosh(y_ab));
389  } else { // weights are zero
390  //pab = PseudoJet(0.0,0.0,0.0,0.0);
391  pab.reset(0.0, 0.0, 0.0, 0.0);
392  }
393 }
394 
395 
397  switch(_recomb_scheme) {
398  case E_scheme:
399  case BIpt_scheme:
400  case BIpt2_scheme:
401  case WTA_pt_scheme:
402  //case WTA_E_scheme:
403  case WTA_modp_scheme:
404  break;
405  case pt_scheme:
406  case pt2_scheme:
407  {
408  // these schemes (as in the ktjet implementation) need massless
409  // initial 4-vectors with essentially E=|p|.
410  double newE = sqrt(p.perp2()+p.pz()*p.pz());
411  p.reset_momentum(p.px(), p.py(), p.pz(), newE);
412  // FJ2.x version
413  // int user_index = p.user_index();
414  // p = PseudoJet(p.px(), p.py(), p.pz(), newE);
415  // p.set_user_index(user_index);
416  }
417  break;
418  case Et_scheme:
419  case Et2_scheme:
420  {
421  // these schemes (as in the ktjet implementation) need massless
422  // initial 4-vectors with essentially E=|p|.
423  double rescale = p.E()/sqrt(p.perp2()+p.pz()*p.pz());
424  p.reset_momentum(rescale*p.px(), rescale*p.py(), rescale*p.pz(), p.E());
425  // FJ2.x version
426  // int user_index = p.user_index();
427  // p = PseudoJet(rescale*p.px(), rescale*p.py(), rescale*p.pz(), p.E());
428  // p.set_user_index(user_index);
429  }
430  break;
431  default:
432  ostringstream err;
433  err << "DefaultRecombiner: unrecognized recombination scheme "
434  << _recomb_scheme;
435  throw Error(err.str());
436  }
437 }
438 
440  throw Error("set_ghost_separation_scale not supported");
441 }
442 
443 
444 
445 //-------------------------------------------------------------------------------
446 // helper functions to build a jet made of pieces
447 //
448 // This is the extended version with support for a user-defined
449 // recombination-scheme
450 // -------------------------------------------------------------------------------
451 
452 // build a "CompositeJet" from the vector of its pieces
453 //
454 // the user passes the reciombination scheme used to "sum" the pieces.
455 PseudoJet join(const vector<PseudoJet> & pieces, const JetDefinition::Recombiner & recombiner){
456  // compute the total momentum
457  //--------------------------------------------------
458  PseudoJet result; // automatically initialised to 0
459  if (pieces.size()>0){
460  result = pieces[0];
461  for (unsigned int i=1; i<pieces.size(); i++)
462  recombiner.plus_equal(result, pieces[i]);
463  }
464 
465  // attach a CompositeJetStructure to the result
466  //--------------------------------------------------
467  CompositeJetStructure *cj_struct = new CompositeJetStructure(pieces, &recombiner);
468 
470 
471  return result;
472 }
473 
474 // build a "CompositeJet" from a single PseudoJet
475 PseudoJet join(const PseudoJet & j1,
476  const JetDefinition::Recombiner & recombiner){
477  return join(vector<PseudoJet>(1,j1), recombiner);
478 }
479 
480 // build a "CompositeJet" from two PseudoJet
481 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2,
482  const JetDefinition::Recombiner & recombiner){
483  vector<PseudoJet> pieces;
484  pieces.push_back(j1);
485  pieces.push_back(j2);
486  return join(pieces, recombiner);
487 }
488 
489 // build a "CompositeJet" from 3 PseudoJet
490 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3,
491  const JetDefinition::Recombiner & recombiner){
492  vector<PseudoJet> pieces;
493  pieces.push_back(j1);
494  pieces.push_back(j2);
495  pieces.push_back(j3);
496  return join(pieces, recombiner);
497 }
498 
499 // build a "CompositeJet" from 4 PseudoJet
500 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, const PseudoJet & j4,
501  const JetDefinition::Recombiner & recombiner){
502  vector<PseudoJet> pieces;
503  pieces.push_back(j1);
504  pieces.push_back(j2);
505  pieces.push_back(j3);
506  pieces.push_back(j4);
507  return join(pieces, recombiner);
508 }
509 
510 
511 
512 
513 FASTJET_END_NAMESPACE
pt weighted recombination of y,phi (and summing of pt&#39;s), with no preprocessing
void set_recombination_scheme(RecombinationScheme)
set the recombination scheme to the one provided
a version of cambridge with a special distance measure for particles whose pt is < extra_param(); thi...
like the k_t but with distance measures dij = min(kti^{2p},ktj^{2p}) Delta R_{ij}^2 / R^2 diB = 1/kti...
virtual void preprocess(PseudoJet &p) const override
routine called to preprocess each input jet (to make all input jets compatible with the scheme requir...
bool is_spherical() const
returns true if the jet definition involves an algorithm intended for use on a spherical geometry (e...
void delete_recombiner_when_unused()
calling this tells the JetDefinition to handle the deletion of the recombiner when it is no longer us...
summing the 4-momenta
double m2() const
returns the squared invariant mass // like CLHEP
Definition: PseudoJet.hh:150
const Plugin * plugin() const
return a pointer to the plugin
JetAlgorithm jet_algorithm() const
return information about the definition...
virtual void recombine(const PseudoJet &pa, const PseudoJet &pb, PseudoJet &pab) const override
recombine pa and pb and put result into pab
the value for the jet algorithm in a JetDefinition for which no algorithm has yet been defined ...
pt^2 weighted recombination of y,phi (and summing of pt&#39;s) no preprocessing
void delete_plugin_when_unused()
calling this causes the JetDefinition to handle the deletion of the plugin when it is no longer used ...
pt^2 weighted recombination of y,phi (and summing of pt&#39;s) with preprocessing to make things massless...
the e+e- genkt algorithm (R > 2 and p=1 gives ee_kt)
std::string description_no_recombiner() const
returns a description not including the recombiner information
void reset_PtYPhiM(double pt_in, double y_in, double phi_in, double m_in=0.0)
reset the PseudoJet according to the specified pt, rapidity, azimuth and mass (also resetting indices...
Definition: PseudoJet.hh:280
the longitudinally invariant kt algorithm
static std::string algorithm_description(const JetAlgorithm jet_alg)
a short textual description of the algorithm jet_alg
pt-based Winner-Takes-All (WTA) recombination: the result of the recombination has the rapidity...
The structure for a jet made of pieces.
double pt() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:141
void set_extra_param(double xtra_param)
(re)set the general purpose extra parameter
double modp2() const
return the squared 3-vector modulus = px^2+py^2+pz^2
Definition: PseudoJet.hh:165
any plugin algorithm supplied by the user
mod-p-based Winner-Takes-All (WTA) recombination: the result of the recombination gets the 3-vector d...
the plugin has been used...
A class that will provide the recombination scheme facilities and/or allow a user to extend these fac...
double perp() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:145
const Recombiner * recombiner() const
returns a pointer to the currently defined recombiner.
pt^2 weighted recombination of y,phi (and summing of pt&#39;s) with preprocessing to make things massless...
virtual std::string description() const =0
return a textual description of the jet-definition implemented in this plugin
base class corresponding to errors that can be thrown by FastJet
Definition: Error.hh:47
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP) ...
Definition: PseudoJet.hh:951
RecombinationScheme
The various recombination schemes.
double perp2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:143
static const double max_allowable_R
R values larger than max_allowable_R are not allowed.
an implementation of C++0x shared pointers (or boost&#39;s)
Definition: SharedPtr.hh:121
pt weighted recombination of y,phi (and summing of pt&#39;s) with preprocessing to make things massless b...
double phi() const
returns phi (in the range 0..2pi)
Definition: PseudoJet.hh:110
double rap() const
returns the rapidity or some large value when the rapidity is infinite
Definition: PseudoJet.hh:125
virtual std::string description() const =0
return a textual description of the recombination scheme implemented here
pt weighted recombination of y,phi (and summing of pt&#39;s) with preprocessing to make things massless b...
double modp() const
return the 3-vector modulus = sqrt(px^2+py^2+pz^2)
Definition: PseudoJet.hh:167
static unsigned int n_parameters_for_algorithm(const JetAlgorithm jet_alg)
the number of parameters associated to a given jet algorithm
void reset(double px, double py, double pz, double E)
reset the 4-momentum according to the supplied components and put the user and history indices back t...
Definition: PseudoJet.hh:957
virtual bool is_spherical() const
returns true if the plugin implements an algorithm intended for use on a spherical geometry (e...
void set_recombiner(const Recombiner *recomb)
set the recombiner class to the one provided
like the k_t but with distance measures dij = min(1/kti^2,1/ktj^2) Delta R_{ij}^2 / R^2 diB = 1/kti^2...
double pt2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:139
the e+e- kt algorithm
virtual void set_ghost_separation_scale(double scale) const
set the ghost separation scale for passive area determinations in future runs (strictly speaking that...
void set_structure_shared_ptr(const SharedPtr< PseudoJetStructureBase > &structure)
set the associated structure
Definition: PseudoJet.cc:468
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
bool has_same_recombiner(const JetDefinition &other_jd) const
returns true if the current jet definitions shares the same recombiner as the one passed as an argume...
std::string description() const
return a textual description of the current jet definition
An abstract base class that will provide the recombination scheme facilities and/or allow a user to e...
JetAlgorithm
the various families of jet-clustering algorithm
Strategy
the various options for the algorithmic strategy to adopt in clustering events with kt and cambridge ...
virtual std::string description() const override
return a textual description of the recombination scheme implemented here
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
void plus_equal(PseudoJet &pa, const PseudoJet &pb) const
pa += pb in the given recombination scheme.
void reset_momentum(double px, double py, double pz, double E)
reset the 4-momentum according to the supplied components but leave all other information (indices...
Definition: PseudoJet.hh:966
for the user&#39;s external scheme
class that is intended to hold a full definition of the jet clusterer