29 #include "fastjet/JetDefinition.hh"
30 #include "fastjet/Error.hh"
31 #include "fastjet/CompositeJetStructure.hh"
34 FASTJET_BEGIN_NAMESPACE
38 const double JetDefinition::max_allowable_R = 1000.0;
48 _jet_algorithm(jet_algorithm_in), _Rparam(R_in), _strategy(strategy_in) {
63 oss <<
"Requested R = " << R_in <<
" for jet definition is larger than max_allowable_R = " <<
max_allowable_R;
64 throw Error(oss.str());
70 switch (jet_algorithm_in) {
72 if (nparameters != 0) {
74 oss <<
"ee_kt_algorithm should be constructed with 0 parameters but was called with "
75 << nparameters <<
" parameter(s)\n";
76 throw Error(oss.str());
81 if (nparameters != 2) {
83 oss <<
"(ee_)genkt_algorithm should be constructed with 2 parameters but was called with "
84 << nparameters <<
" parameter(s)\n";
85 throw Error(oss.str());
89 if (nparameters != 1) {
91 oss <<
"The jet algorithm you requested ("
92 << jet_algorithm_in <<
") should be constructed with 1 parameter but was called with "
93 << nparameters <<
" parameter(s)\n";
94 throw Error(oss.str());
113 name <<
"Longitudinally invariant kt algorithm with R = " << R();
116 name <<
"Longitudinally invariant Cambridge/Aachen algorithm with R = "
120 name <<
"Longitudinally invariant anti-kt algorithm with R = "
124 name <<
"Longitudinally invariant generalised kt algorithm with R = "
125 << R() <<
", p = " << extra_param();
128 name <<
"Longitudinally invariant Cambridge/Aachen algorithm with R = "
129 << R() <<
"and a special hack whereby particles with kt < "
130 << extra_param() <<
"are treated as passive ghosts";
132 name <<
"e+e- kt (Durham) algorithm (NB: no R)";
135 name <<
"e+e- generalised kt algorithm with R = "
136 << R() <<
", p = " << extra_param();
139 name <<
"uninitialised JetDefinition (jet_algorithm=undefined_jet_algorithm)" ;
141 throw Error(
"JetDefinition::description(): unrecognized jet_algorithm");
152 if (_recombiner_shared()) _recombiner_shared.reset();
163 if (other_jd.recombination_scheme() != scheme)
return false;
174 if (_recombiner == 0){
175 throw Error(
"tried to call JetDefinition::delete_recombiner_when_unused() for a JetDefinition without a user-defined recombination scheme");
178 _recombiner_shared.reset(_recombiner);
185 throw Error(
"tried to call JetDefinition::delete_plugin_when_unused() for a JetDefinition without a plugin");
188 _plugin_shared.reset(_plugin);
194 switch(_recomb_scheme) {
196 return "E scheme recombination";
198 return "pt scheme recombination";
200 return "pt2 scheme recombination";
202 return "Et scheme recombination";
204 return "Et2 scheme recombination";
206 return "boost-invariant pt scheme recombination";
208 return "boost-invariant pt2 scheme recombination";
211 err <<
"DefaultRecombiner: unrecognized recombination scheme "
213 throw Error(err.str());
222 double weighta, weightb;
224 switch(_recomb_scheme) {
229 pab.
reset(pa.px()+pb.px(),
245 weighta = pa.
perp2();
246 weightb = pb.
perp2();
250 err <<
"DefaultRecombiner: unrecognized recombination scheme "
252 throw Error(err.str());
255 double perp_ab = pa.
perp() + pb.
perp();
256 if (perp_ab != 0.0) {
257 double y_ab = (weighta * pa.
rap() + weightb * pb.
rap())/(weighta+weightb);
260 double phi_a = pa.
phi(), phi_b = pb.
phi();
261 if (phi_a - phi_b > pi) phi_b += twopi;
262 if (phi_a - phi_b < -pi) phi_b -= twopi;
263 double phi_ab = (weighta * phi_a + weightb * phi_b)/(weighta+weightb);
273 pab.
reset(0.0, 0.0, 0.0, 0.0);
279 switch(_recomb_scheme) {
289 double newE = sqrt(p.
perp2()+p.pz()*p.pz());
302 double rescale = p.E()/sqrt(p.
perp2()+p.pz()*p.pz());
303 p.
reset_momentum(rescale*p.px(), rescale*p.py(), rescale*p.pz(), p.E());
312 err <<
"DefaultRecombiner: unrecognized recombination scheme "
314 throw Error(err.str());
319 throw Error(
"set_ghost_separation_scale not supported");
338 if (pieces.size()>0){
340 for (
unsigned int i=1; i<pieces.size(); i++)
356 return join(vector<PseudoJet>(1,j1),
recombiner);
362 vector<PseudoJet> pieces;
363 pieces.push_back(j1);
364 pieces.push_back(j2);
371 vector<PseudoJet> pieces;
372 pieces.push_back(j1);
373 pieces.push_back(j2);
374 pieces.push_back(j3);
381 vector<PseudoJet> pieces;
382 pieces.push_back(j1);
383 pieces.push_back(j2);
384 pieces.push_back(j3);
385 pieces.push_back(j4);
392 FASTJET_END_NAMESPACE
const Recombiner * recombiner() const
return a pointer to the currently defined recombiner.
@ BIpt_scheme
pt weighted recombination of y,phi (and summing of pt's), with no preprocessing
@ BIpt2_scheme
pt^2 weighted recombination of y,phi (and summing of pt's) no preprocessing
virtual std::string description() const =0
return a textual description of the recombination scheme implemented here
@ Et_scheme
pt weighted recombination of y,phi (and summing of pt's) with preprocessing to make things massless b...
static const double max_allowable_R
R values larger than max_allowable_R are not allowed.
@ E_scheme
summing the 4-momenta
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...
void set_recombination_scheme(RecombinationScheme)
set the recombination scheme to the one provided
virtual void recombine(const PseudoJet &pa, const PseudoJet &pb, PseudoJet &pab) const
recombine pa and pb and put result into pab
@ cambridge_algorithm
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
virtual void preprocess(PseudoJet &p) const
routine called to preprocess each input jet (to make all input jets compatible with the scheme requir...
double perp() const
returns the scalar transverse momentum
std::string description() const
return a textual description of the current jet definition
@ pt2_scheme
pt^2 weighted recombination of y,phi (and summing of pt's) with preprocessing to make things massless...
void plus_equal(PseudoJet &pa, const PseudoJet &pb) const
pa += pb in the given recombination scheme.
JetAlgorithm jet_algorithm() const
return information about the definition...
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
@ genkt_algorithm
like the k_t but with distance measures dij = min(kti^{2p},ktj^{2p}) Delta R_{ij}^2 / R^2 diB = 1/kti...
bool has_same_recombiner(const JetDefinition &other_jd) const
returns true if the current jet definitions shares the same recombiner as teh one passed as an argume...
void set_extra_param(double xtra_param)
(re)set the general purpose extra parameter
@ antikt_algorithm
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
@ ee_genkt_algorithm
the e+e- genkt algorithm (R > 2 and p=1 gives ee_kt)
@ external_scheme
for the user's external scheme
double rap() const
returns the rapidity or some large value when the rapidity is infinite
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,...
virtual std::string description() const =0
return a textual description of the jet-definition implemented in this plugin
Strategy
the various options for the algorithmic strategy to adopt in clustering events with kt and cambridge ...
@ Et2_scheme
pt^2 weighted recombination of y,phi (and summing of pt's) with preprocessing to make things massless...
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...
double phi() const
returns phi (in the range 0..2pi)
RecombinationScheme
the various recombination schemes
@ undefined_jet_algorithm
the value for the jet algorithm in a JetDefinition for which no algorithm has yet been defined
@ plugin_strategy
the plugin has been used...
void delete_recombiner_when_unused()
calling this tells the JetDefinition to handle the deletion of the recombiner when it is no longer us...
double perp2() const
returns the squared transverse momentum
@ plugin_algorithm
any plugin algorithm supplied by the user
virtual std::string description() const
return a textual description of the recombination scheme implemented here
const Plugin * plugin() const
return a pointer to the plugin
@ pt_scheme
pt weighted recombination of y,phi (and summing of pt's) with preprocessing to make things massless b...
@ ee_kt_algorithm
the e+e- kt algorithm
@ cambridge_for_passive_algorithm
a version of cambridge with a special distance measure for particles whose pt is < extra_param()
@ kt_algorithm
the longitudinally invariant kt algorithm
void delete_plugin_when_unused()
allows to let the JetDefinition handle the deletion of the plugin when it is no longer used