FastJet  3.0.6
fastjet_timing_plugins.cc
1 //STARTHEADER
2 // $Id: fastjet_timing_plugins.cc 3209 2013-09-15 20:03:30Z salam $
3 //
4 // Copyright (c) 2005-2011, 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 and are described in hep-ph/0512210. If you use
16 // FastJet as part of work towards a scientific publication, please
17 // include a citation to the FastJet paper.
18 //
19 // FastJet is distributed in the hope that it will be useful,
20 // but WITHOUT ANY WARRANTY; without even the implied warranty of
21 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 // GNU General Public License for more details.
23 //
24 // You should have received a copy of the GNU General Public License
25 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
26 //----------------------------------------------------------------------
27 //ENDHEADER
28 
29 
30 //----------------------------------------------------------------------
31 /// fastjet_timing.cc: Program to help time and test the fastjet package
32 ///
33 /// It reads files containing multiple events in the format
34 /// p1x p1y p1z E1
35 /// p2x p2y p2z E2
36 /// ...
37 /// #END
38 ///
39 /// An example input file containing 10 events is included as
40 /// data/Pythia-PtMin1000-LHC-10ev.dat
41 ///
42 /// Usage:
43 /// fastjet_timing [-strategy NUMBER] [-repeat nrepeats] [-massive] \
44 /// [-combine nevents] [-r Rparameter] [-incl ptmin] [...] \
45 /// < data_file
46 ///
47 /// where the clustering can be repeated to aid timing and multiple
48 /// events can be combined to get to larger multiplicities. Some options:
49 ///
50 /// Options for reading
51 /// -------------------
52 ///
53 /// -nev n number of events to run
54 ///
55 /// -combine n for combining multiple events from the data file in order
56 /// to get a single high-multipicity event to run.
57 ///
58 /// -massless read in only the 3-momenta and deduce energies assuming
59 /// that particles are massless
60 ///
61 /// -dense adds dense ghost coverage
62 ///
63 /// -repeat n repeats each event n times
64 ///
65 /// -nhardest n keep only the n hardest particles in the event
66 ///
67 /// -file name read from the corresponding file rather than stdin.
68 /// (The file will be reopened for each new jet alg.; in
69 /// constrast, if you use stdin, each new alg will take a
70 /// new event).
71 ///
72 /// Output Options
73 /// --------------
74 ///
75 /// -incl ptmin output of all inclusive jets with pt > ptmin is obtained
76 /// with the -incl option.
77 ///
78 /// -repeat-incl ptmin
79 /// same as -incl ptmin but do it for each repetition
80 /// of the clustering
81 ///
82 /// -excld dcut output of all exclusive jets as obtained in a clustering
83 /// with dcut
84 ///
85 /// -excly ycut output of all exclusive jets as obtained in a clustering
86 /// with ycut
87 ///
88 /// -excln n output of clustering to n exclusive jets
89 ///
90 /// -ee-print print things as px,py,pz,E
91 ///
92 /// -get-all-dij print out all dij values
93 /// -get-all-yij print out all yij values
94 ///
95 /// -const show jet constituents (works with excl jets)
96 ///
97 /// -write for writing out detailed clustering sequence (valuable
98 /// for testing purposes)
99 ///
100 /// -unique_write writes out the sequence of dij's according to the
101 /// "unique_history_order" (useful for verifying consistency
102 /// between different clustering strategies).
103 ///
104 /// -root file sends output to file that can be read in with the script in
105 /// root/ so as to show a lego-plot of the event
106 ///
107 /// -cones show extra info about internal steps for SISCone
108 ///
109 /// -area calculate areas. Additional options include
110 /// -area:active
111 /// -area:passive
112 /// -area:explicit
113 /// -area:voronoi Rfact
114 /// -area:repeat nrepeat
115 /// -ghost-area area
116 /// -ghost-maxrap maxrap
117 /// -area:fj2 place ghosts as in fj2
118 ///
119 /// -bkgd calculate the background density. Additional options include
120 /// -bkgd:csab use the old ClusterSequenceAreaBase methods
121 /// -bkgd:jetmedian use the new JetMedianBackgroundEstimator class
122 /// -bkgd:fj2 force jetmedian to calculate sigma as in fj2
123 /// -bkgd:gridmedian use GridMedianBackgroundEstimator with grid up to ghost_maxrap-ktR and grid spacing of 2ktR
124 ///
125 /// Algorithms
126 /// ----------
127 /// -all-algs runs all algorithms
128 ///
129 /// -kt switch to the longitudinally invariant kt algorithm
130 /// Note: this is the default one.
131 ///
132 /// -cam switch to the inclusive Cambridge/Aachen algorithm --
133 /// note that the option -excld dcut provides a clustering
134 /// up to the dcut which is the minimum squared
135 /// distance between any pair of jets.
136 ///
137 /// -antikt switch to the anti-kt clustering algorithm
138 ///
139 /// -genkt switch to the genkt algorithm
140 /// you can provide the parameter of the alg as an argument to
141 /// -genkt (1 by default)
142 ///
143 /// -eekt switch to the e+e- kt algorithm
144 ///
145 /// -eegenkt switch to the genkt algorithm
146 /// you can provide the parameter of the alg as an argument to
147 /// -ee_genkt (1 by default)
148 ///
149 /// plugins (don't delete this line)
150 ///
151 /// -pxcone switch to the PxCone jet algorithm
152 ///
153 /// -siscone switch to the SISCone jet algorithm (seedless cones)
154 /// -sisconespheri switch to the Spherical SISCone jet algorithm (seedless cones)
155 ///
156 /// -midpoint switch to CDF's midpoint code
157 /// -jetclu switch to CDF's jetclu code
158 ///
159 /// -d0runipre96cone switch to the D0RunIpre96Cone plugin
160 /// -d0runicone switch to the D0RunICone plugin
161 ///
162 /// -d0runiicone switch to D0's run II midpoint cone
163 ///
164 /// -trackjet switch to the TrackJet plugin
165 ///
166 /// -atlascone switch to the ATLASCone plugin
167 ///
168 /// -eecambridge switch to the EECambridge plugin
169 ///
170 /// -jade switch to the Jade plugin
171 ///
172 /// -cmsiterativecone switch to the CMSIterativeCone plugin
173 ///
174 /// -gridjet switch to the GridJet plugin
175 ///
176 /// end of plugins (don't delete this line)
177 ///
178 ///
179 /// Options for running algs
180 /// ------------------------
181 ///
182 /// -r sets the radius of the jet algorithm (default = 1.0)
183 ///
184 /// -overlap | -f sets the overlap fraction in cone algs with split-merge
185 ///
186 /// -seed sets the seed threshold
187 ///
188 /// -strategy N indicate stratgey from the enum fastjet::Strategy (see
189 /// fastjet/JetDefinition.hh).
190 ///
191 
192 
193 #include "fastjet/ClusterSequenceArea.hh"
194 #include "fastjet/tools/JetMedianBackgroundEstimator.hh"
195 #include "fastjet/tools/GridMedianBackgroundEstimator.hh"
196 #include "fastjet/Selector.hh"
197 #include<iostream>
198 #include<sstream>
199 #include<fstream>
200 #include<valarray>
201 #include<vector>
202 #include <cstdlib>
203 //#include<cstddef> // for size_t
204 #include "CmdLine.hh"
205 
206 // get info on how fastjet was configured
207 #include "fastjet/config.h"
208 
209 // include the installed plugins (don't delete this line)
210 #ifdef FASTJET_ENABLE_PLUGIN_SISCONE
211 #include "fastjet/SISConePlugin.hh"
212 #include "fastjet/SISConeSphericalPlugin.hh"
213 #endif
214 #ifdef FASTJET_ENABLE_PLUGIN_CDFCONES
215 #include "fastjet/CDFMidPointPlugin.hh"
216 #include "fastjet/CDFJetCluPlugin.hh"
217 #endif
218 #ifdef FASTJET_ENABLE_PLUGIN_PXCONE
219 #include "fastjet/PxConePlugin.hh"
220 #endif
221 #ifdef FASTJET_ENABLE_PLUGIN_D0RUNIICONE
222 #include "fastjet/D0RunIIConePlugin.hh"
223 #endif
224 #ifdef FASTJET_ENABLE_PLUGIN_TRACKJET
225 #include "fastjet/TrackJetPlugin.hh"
226 #endif
227 #ifdef FASTJET_ENABLE_PLUGIN_ATLASCONE
228 #include "fastjet/ATLASConePlugin.hh"
229 #endif
230 #ifdef FASTJET_ENABLE_PLUGIN_EECAMBRIDGE
231 #include "fastjet/EECambridgePlugin.hh"
232 #endif
233 #ifdef FASTJET_ENABLE_PLUGIN_JADE
234 #include "fastjet/JadePlugin.hh"
235 #endif
236 #ifdef FASTJET_ENABLE_PLUGIN_CMSITERATIVECONE
237 #include "fastjet/CMSIterativeConePlugin.hh"
238 #endif
239 #ifdef FASTJET_ENABLE_PLUGIN_D0RUNICONE
240 #include "fastjet/D0RunIpre96ConePlugin.hh"
241 #include "fastjet/D0RunIConePlugin.hh"
242 #endif
243 #ifdef FASTJET_ENABLE_PLUGIN_GRIDJET
244 #include "fastjet/GridJetPlugin.hh"
245 #endif
246 // end of installed plugins inclusion (don't delete this line)
247 
248 using namespace std;
249 
250 // to avoid excessive typing, define an abbreviation for the
251 // fastjet namespace
252 namespace fj = fastjet;
253 
254 inline double pow2(const double x) {return x*x;}
255 
256 // pretty print the jets and their subjets
257 void print_jets_and_sub (const vector<fj::PseudoJet> & jets, double dcut);
258 
259 string rootfile;
260 CmdLine * cmdline_p;
261 
262 bool do_areas;
263 
264 /// sort and pretty print jets, with exact behaviour depending on
265 /// whether ee_print is true or not
266 bool ee_print = false;
267 void print_jets(const vector<fj::PseudoJet> & jets, bool show_const = false);
268 
269 bool found_unavailable = false;
270 void is_unavailable(const string & algname) {
271  cerr << algname << " requested, but not available for this compilation" << endl;
272  found_unavailable = true;
273  //exit(0);
274 }
275 
276 
277 /// a program to test and time a range of algorithms as implemented or
278 /// wrapped in fastjet
279 int main (int argc, char ** argv) {
280 
282 
283  CmdLine cmdline(argc,argv);
284  cmdline_p = &cmdline;
285  // allow the use to specify the fj::Strategy either through the
286  // -clever or the -strategy options (both will take numerical
287  // values); the latter will override the former.
288  fj::Strategy strategy = fj::Strategy(cmdline.int_val("-strategy",
289  cmdline.int_val("-clever", fj::Best)));
290  int repeat = cmdline.int_val("-repeat",1);
291  int combine = cmdline.int_val("-combine",1);
292  bool write = cmdline.present("-write");
293  bool unique_write = cmdline.present("-unique_write");
294  bool hydjet = cmdline.present("-hydjet");
295  double ktR = cmdline.double_val("-r",1.0);
296  ktR = cmdline.double_val("-R",ktR); // allow -r and -R
297  double inclkt = cmdline.double_val("-incl",-1.0);
298  double repeatinclkt = cmdline.double_val("-repeat-incl",-1.0);
299  int excln = cmdline.int_val ("-excln",-1);
300  double excld = cmdline.double_val("-excld",-1.0);
301  double excly = cmdline.double_val("-excly",-1.0);
302  ee_print = cmdline.present("-ee-print");
303  bool get_all_dij = cmdline.present("-get-all-dij");
304  bool get_all_yij = cmdline.present("-get-all-yij");
305  double subdcut = cmdline.double_val("-subdcut",-1.0);
306  double rapmax = cmdline.double_val("-rapmax",1.0e305);
307  if (cmdline.present("-etamax")) {
308  cerr << "WARNING: -etamax options actually sets maximum rapidity (and overrides -rapmax)\n";
309  rapmax = cmdline.double_val("-etamax");
310  }
311  bool show_constituents = cmdline.present("-const");
312  bool massless = cmdline.present("-massless");
313  int nev = cmdline.int_val("-nev",1);
314  bool add_dense_coverage = cmdline.present("-dense");
315  double ghost_maxrap = cmdline.value("-ghost-maxrap",5.0);
316  bool all_algs = cmdline.present("-all-algs");
317 
318  fj::Selector particles_sel = (cmdline.present("-nhardest"))
319  ? fj::SelectorNHardest(cmdline.value<unsigned int>("-nhardest"))
320  : fj::SelectorIdentity();
321 
322  do_areas = cmdline.present("-area");
323  fj::AreaDefinition area_def;
324  if (do_areas) {
325  assert(!write); // it's incompatible
326  fj::GhostedAreaSpec ghost_spec(ghost_maxrap,
327  cmdline.value("-area:repeat", 1),
328  cmdline.value("-ghost-area", 0.01));
329  if (cmdline.present("-area:fj2")) ghost_spec.set_fj2_placement(true);
330  if (cmdline.present("-area:explicit")) {
331  area_def = fj::AreaDefinition(fj::active_area_explicit_ghosts, ghost_spec);
332  } else if (cmdline.present("-area:passive")) {
333  area_def = fj::AreaDefinition(fj::passive_area, ghost_spec);
334  } else if (cmdline.present("-area:voronoi")) {
335  double Rfact = cmdline.value<double>("-area:voronoi");
336  area_def = fj::AreaDefinition(fj::voronoi_area,
337  fj::VoronoiAreaSpec(Rfact));
338  } else {
339  cmdline.present("-area:active"); // allow, but do not require, arg
340  area_def = fj::AreaDefinition(fj::active_area, ghost_spec);
341  }
342  }
343  bool do_bkgd = cmdline.present("-bkgd"); // background estimation
344  bool do_bkgd_csab = false, do_bkgd_jetmedian = false, do_bkgd_fj2 = false;
345  bool do_bkgd_gridmedian = false;
346  fj::Selector bkgd_range;
347  if (do_bkgd) {
348  bkgd_range = fj::SelectorAbsRapMax(ghost_maxrap - ktR);
349  if (cmdline.present("-bkgd:csab")) {do_bkgd_csab = true;}
350  else if (cmdline.present("-bkgd:jetmedian")) {do_bkgd_jetmedian = true;
351  do_bkgd_fj2 = cmdline.present("-bkgd:fj2");
352  } else if (cmdline.present("-bkgd:gridmedian")) {do_bkgd_gridmedian = true;
353  } else {
354  throw fj::Error("with the -bkgd option, some particular background must be specified (csab or jetmedian)");
355  }
356  assert(do_areas || do_bkgd_gridmedian);
357  }
358 
359  bool show_cones = cmdline.present("-cones"); // only works for siscone
360 
361  // for cone algorithms
362  // allow -f and -overlap
363  double overlap_threshold = cmdline.double_val("-overlap",0.5);
364  overlap_threshold = cmdline.double_val("-f",overlap_threshold);
365  double seed_threshold = cmdline.double_val("-seed",1.0);
366 
367  // for ee algorithms, allow to specify ycut
368  double ycut = cmdline.double_val("-ycut",0.08);
369 
370  // for printing jets to a file for reading by root
371  rootfile = cmdline.value<string>("-root","");
372 
373  // out default scheme is the E_scheme
375 
376  // The following option causes the Cambridge algo to be used.
377  // Note that currently the only output that works sensibly here is
378  // "-incl 0"
379  vector<fj::JetDefinition> jet_defs;
380  if (all_algs || cmdline.present("-cam") || cmdline.present("-CA")) {
381  jet_defs.push_back( fj::JetDefinition(fj::cambridge_algorithm, ktR, scheme, strategy));
382  }
383  if (all_algs || cmdline.present("-antikt")) {
384  jet_defs.push_back( fj::JetDefinition(fj::antikt_algorithm, ktR, scheme, strategy));
385  }
386  if (all_algs || cmdline.present("-genkt")) {
387  double p;
388  if (cmdline.present("-genkt")) p = cmdline.value<double>("-genkt");
389  else p = -0.5;
390  jet_defs.push_back( fj::JetDefinition(fj::genkt_algorithm, ktR, p, scheme, strategy));
391  }
392  if (all_algs || cmdline.present("-eekt")) {
393  jet_defs.push_back( fj::JetDefinition(fj::ee_kt_algorithm));
394  }
395  if (all_algs || cmdline.present("-eegenkt")) {
396  double p;
397  if (cmdline.present("-eegenkt")) p = cmdline.value<double>("-eegenkt");
398  else p = -0.5;
399  jet_defs.push_back( fj::JetDefinition(fj::ee_genkt_algorithm, ktR, p, scheme, strategy));
400 
401 // checking if one asks to run a plugin (don't delete this line)
402  }
403  if (all_algs || cmdline.present("-midpoint")) {
404 #ifdef FASTJET_ENABLE_PLUGIN_CDFCONES
405  typedef fj::CDFMidPointPlugin MPPlug; // for brevity
406  double cone_area_fraction = 1.0;
407  int max_pair_size = 2;
408  int max_iterations = 100;
409  MPPlug::SplitMergeScale sm_scale = MPPlug::SM_pt;
410  if (cmdline.present("-sm-pttilde")) sm_scale = MPPlug::SM_pttilde;
411  if (cmdline.present("-sm-pt")) sm_scale = MPPlug::SM_pt; // default
412  if (cmdline.present("-sm-mt")) sm_scale = MPPlug::SM_mt;
413  if (cmdline.present("-sm-Et")) sm_scale = MPPlug::SM_Et;
414  jet_defs.push_back( fj::JetDefinition( new fj::CDFMidPointPlugin (
415  seed_threshold, ktR,
416  cone_area_fraction, max_pair_size,
417  max_iterations, overlap_threshold,
418  sm_scale)));
419 #else // FASTJET_ENABLE_PLUGIN_CDFCONES
420  is_unavailable("MidPoint");
421 #endif // FASTJET_ENABLE_PLUGIN_CDFCONES
422  }
423  if (all_algs || cmdline.present("-pxcone")) {
424 #ifdef FASTJET_ENABLE_PLUGIN_PXCONE
425  double min_jet_energy = 5.0;
426  jet_defs.push_back( fj::JetDefinition( new fj::PxConePlugin (
427  ktR, min_jet_energy,
428  overlap_threshold)));
429 #else // FASTJET_ENABLE_PLUGIN_PXCONE
430  is_unavailable("PxCone");
431 #endif // FASTJET_ENABLE_PLUGIN_PXCONE
432  }
433  if (all_algs || cmdline.present("-jetclu")) {
434 #ifdef FASTJET_ENABLE_PLUGIN_CDFCONES
435  jet_defs.push_back( fj::JetDefinition( new fj::CDFJetCluPlugin (
436  ktR, overlap_threshold, seed_threshold)));
437 #else // FASTJET_ENABLE_PLUGIN_CDFCONES
438  is_unavailable("JetClu");
439 #endif // FASTJET_ENABLE_PLUGIN_CDFCONES
440  }
441  if (all_algs || cmdline.present("-siscone") || cmdline.present("-sisconespheri")) {
442 #ifdef FASTJET_ENABLE_PLUGIN_SISCONE
443  typedef fj::SISConePlugin SISPlug; // for brevity
444  int npass = cmdline.value("-npass",0);
445  if (all_algs || cmdline.present("-siscone")) {
446  double sisptmin = cmdline.value("-sisptmin",0.0);
447  SISPlug * plugin = new SISPlug (ktR, overlap_threshold,npass,sisptmin);
448  if (cmdline.present("-sm-pt")) plugin->set_split_merge_scale(SISPlug::SM_pt);
449  if (cmdline.present("-sm-mt")) plugin->set_split_merge_scale(SISPlug::SM_mt);
450  if (cmdline.present("-sm-Et")) plugin->set_split_merge_scale(SISPlug::SM_Et);
451  if (cmdline.present("-sm-pttilde")) plugin->set_split_merge_scale(SISPlug::SM_pttilde);
452  // cause it to use the jet-definition's own recombiner
453  plugin->set_use_jet_def_recombiner(true);
454  jet_defs.push_back( fj::JetDefinition(plugin));
455  }
456  if (all_algs || cmdline.present("-sisconespheri")) {
457  double sisEmin = cmdline.value("-sisEmin",0.0);
458  fj::SISConeSphericalPlugin * plugin =
459  new fj::SISConeSphericalPlugin(ktR, overlap_threshold,npass,sisEmin);
460  if (cmdline.present("-ghost-sep")) {
461  plugin->set_ghost_separation_scale(cmdline.value<double>("-ghost-sep"));
462  }
463  jet_defs.push_back( fj::JetDefinition(plugin));
464  }
465 #else // FASTJET_ENABLE_PLUGIN_SISCONE
466  is_unavailable("SISCone");
467 #endif // FASTJET_ENABLE_PLUGIN_SISCONE
468  }
469  if (all_algs || cmdline.present("-d0runiicone")) {
470 #ifdef FASTJET_ENABLE_PLUGIN_D0RUNIICONE
471  double min_jet_Et = 6.0; // was 8 GeV in earlier work
472  jet_defs.push_back( fj::JetDefinition(new fj::D0RunIIConePlugin(ktR,min_jet_Et)));
473 #else // FASTJET_ENABLE_PLUGIN_D0RUNIICONE
474  is_unavailable("D0RunIICone");
475 #endif // FASTJET_ENABLE_PLUGIN_D0RUNIICONE
476  }
477  if (all_algs || cmdline.present("-trackjet")) {
478 #ifdef FASTJET_ENABLE_PLUGIN_TRACKJET
479  jet_defs.push_back( fj::JetDefinition(new fj::TrackJetPlugin(ktR)));
480 #else // FASTJET_ENABLE_PLUGIN_TRACKJET
481  is_unavailable("TrackJet");
482 #endif // FASTJET_ENABLE_PLUGIN_TRACKJET
483  }
484  if (all_algs || cmdline.present("-atlascone")) {
485 #ifdef FASTJET_ENABLE_PLUGIN_ATLASCONE
486  jet_defs.push_back( fj::JetDefinition(new fj::ATLASConePlugin(ktR)));
487 #else // FASTJET_ENABLE_PLUGIN_ATLASCONE
488  is_unavailable("ATLASCone");
489 #endif // FASTJET_ENABLE_PLUGIN_ATLASCONE
490  }
491  if (all_algs || cmdline.present("-eecambridge")) {
492 #ifdef FASTJET_ENABLE_PLUGIN_EECAMBRIDGE
493  jet_defs.push_back( fj::JetDefinition(new fj::EECambridgePlugin(ycut)));
494 #else // FASTJET_ENABLE_PLUGIN_EECAMBRIDGE
495  is_unavailable("EECambridge");
496 #endif // FASTJET_ENABLE_PLUGIN_EECAMBRIDGE
497  }
498  if (all_algs || cmdline.present("-jade")) {
499 #ifdef FASTJET_ENABLE_PLUGIN_JADE
500  jet_defs.push_back( fj::JetDefinition(new fj::JadePlugin()));
501 #else // FASTJET_ENABLE_PLUGIN_JADE
502  is_unavailable("Jade");
503 #endif // FASTJET_ENABLE_PLUGIN_JADE
504  }
505  if (all_algs || cmdline.present("-cmsiterativecone")) {
506 #ifdef FASTJET_ENABLE_PLUGIN_CMSITERATIVECONE
507  jet_defs.push_back( fj::JetDefinition(new fj::CMSIterativeConePlugin(ktR,seed_threshold)));
508 #else // FASTJET_ENABLE_PLUGIN_CMSITERATIVECONE
509  is_unavailable("CMSIterativeCone");
510 #endif // FASTJET_ENABLE_PLUGIN_CMSITERATIVECONE
511  }
512  if (all_algs || cmdline.present("-d0runipre96cone")) {
513 #ifdef FASTJET_ENABLE_PLUGIN_D0RUNICONE
514  jet_defs.push_back( fj::JetDefinition(new fj::D0RunIpre96ConePlugin(ktR, seed_threshold, overlap_threshold)));
515 #else // FASTJET_ENABLE_PLUGIN_D0RUNICONE
516  is_unavailable("D0RunIpre96Cone");
517 #endif // FASTJET_ENABLE_PLUGIN_D0RUNICONE
518  }
519  if (all_algs || cmdline.present("-d0runicone")) {
520 #ifdef FASTJET_ENABLE_PLUGIN_D0RUNICONE
521  jet_defs.push_back( fj::JetDefinition(new fj::D0RunIConePlugin(ktR, seed_threshold, overlap_threshold)));
522 #else // FASTJET_ENABLE_PLUGIN_D0RUNICONE
523  is_unavailable("D0RunICone");
524 #endif // FASTJET_ENABLE_PLUGIN_D0RUNICONE
525  }
526  if (all_algs || cmdline.present("-gridjet")) {
527 #ifdef FASTJET_ENABLE_PLUGIN_GRIDJET
528  // we want a grid_ymax of 5.0, but when using R=0.4 (i.e. grid
529  // spacing of 0.8), this leads to 12.5 grid cells; depending on
530  // whether this is 12.499999999999 or 12.5000000....1 this gets
531  // converted either to 12 or 13, making the results sensitive to
532  // rounding errors.
533  //
534  // Instead we therefore take 4.9999999999, which avoids this problem.
535  double grid_ymax = 4.9999999999;
536  jet_defs.push_back( fj::JetDefinition(new fj::GridJetPlugin(grid_ymax, ktR*2.0)));
537 #else // FASTJET_ENABLE_PLUGIN_GRIDJET
538  is_unavailable("GridJet");
539 #endif // FASTJET_ENABLE_PLUGIN_GRIDJET
540 // end of checking if one asks to run a plugin (don't delete this line)
541  }
542  if (all_algs ||
543  cmdline.present("-kt") ||
544  (jet_defs.size() == 0 && !found_unavailable)) {
545  jet_defs.push_back( fj::JetDefinition(fj::kt_algorithm, ktR, strategy));
546  }
547 
548  string filename = cmdline.value<string>("-file", "");
549 
550 
551  if (!cmdline.all_options_used()) {cerr <<
552  "Error: some options were not recognized"<<endl;
553  exit(-1);}
554 
555  for (unsigned idef = 0; idef < jet_defs.size(); idef++) {
556  fj::JetDefinition & jet_def = jet_defs[idef];
557  istream * istr;
558  if (filename == "") istr = &cin;
559  else istr = new ifstream(filename.c_str());
560 
561  for (int iev = 0; iev < nev; iev++) {
562  vector<fj::PseudoJet> jets;
563  vector<fj::PseudoJet> particles;
564  string line;
565  int ndone = 0;
566  while (getline(*istr, line)) {
567  //cout << line<<endl;
568  istringstream linestream(line);
569  if (line == "#END") {
570  ndone += 1;
571  if (ndone == combine) {break;}
572  }
573  if (line.substr(0,1) == "#") {continue;}
574  valarray<double> fourvec(4);
575  if (hydjet) {
576  // special reading from hydjet.txt event record (though actually
577  // this is supposed to be a standard pythia event record, so
578  // being able to read from it is perhaps not so bad an idea...)
579  int ii, istat,id,m1,m2,d1,d2;
580  double mass;
581  linestream >> ii>> istat >> id >> m1 >> m2 >> d1 >> d2
582  >> fourvec[0] >> fourvec[1] >> fourvec[2] >> mass;
583  // current file contains mass of particle as 4th entry
584  if (istat == 1) {
585  fourvec[3] = sqrt(+pow2(fourvec[0])+pow2(fourvec[1])
586  +pow2(fourvec[2])+pow2(mass));
587  }
588  } else {
589  if (massless) {
590  linestream >> fourvec[0] >> fourvec[1] >> fourvec[2];
591  fourvec[3] = sqrt(pow2(fourvec[0])+pow2(fourvec[1])+pow2(fourvec[2]));}
592  else {
593  linestream >> fourvec[0] >> fourvec[1] >> fourvec[2] >> fourvec[3];
594  }
595  }
596  fj::PseudoJet psjet(fourvec);
597  if (abs(psjet.rap()) < rapmax) {particles.push_back(psjet);}
598  }
599 
600  // add a fake underlying event which is very soft, uniformly distributed
601  // in eta,phi so as to allow one to reconstruct the area that is associated
602  // with each jet.
603  if (add_dense_coverage) {
604  fj::GhostedAreaSpec ghosted_area_spec(ghost_maxrap);
605  //fj::GhostedAreaSpec ghosted_area_spec(-2.0,4.0); // asymmetric range
606  // for plots, reduce the scatter default of 1, to avoid "holes"
607  // in the subsequent calorimeter view
608  ghosted_area_spec.set_grid_scatter(0.5);
609  ghosted_area_spec.add_ghosts(particles);
610  //----- old code ------------------
611  // srand(2);
612  // int nphi = 60;
613  // int neta = 100;
614  // double kt = 1e-1;
615  // for (int iphi = 0; iphi<nphi; iphi++) {
616  // for (int ieta = -neta; ieta<neta+1; ieta++) {
617  // double phi = (iphi+0.5) * (fj::twopi/nphi) + rand()*0.001/RAND_MAX;
618  // double eta = ieta * (10.0/neta) + rand()*0.001/RAND_MAX;
619  // kt = 1e-20*(1+rand()*0.1/RAND_MAX);
620  // double pminus = kt*exp(-eta);
621  // double pplus = kt*exp(+eta);
622  // double px = kt*sin(phi);
623  // double py = kt*cos(phi);
624  // //cout << kt<<" "<<eta<<" "<<phi<<"\n";
625  // fj::PseudoJet mom(px,py,0.5*(pplus-pminus),0.5*(pplus+pminus));
626  // particles.push_back(mom);
627  // }
628  // }
629  }
630 
631  // select the particles that pass the selection cut
632  particles = particles_sel(particles);
633 
634  for (int irepeat = 0; irepeat < repeat ; irepeat++) {
635  int nparticles = particles.size();
636  try {
637  auto_ptr<fj::ClusterSequence> clust_seq;
638  if (do_areas) {
639  clust_seq.reset(new fj::ClusterSequenceArea(particles,jet_def,area_def));
640  } else {
641  clust_seq.reset(new fj::ClusterSequence(particles,jet_def,write));
642  }
643 
644  // repetitive output
645  if (repeatinclkt >= 0.0) {
646  vector<fj::PseudoJet> jets_local = sorted_by_pt(clust_seq->inclusive_jets(repeatinclkt));
647  }
648 
649  if (irepeat != 0) {continue;}
650  cout << "iev "<<iev<< ": number of particles = "<< nparticles << endl;
651  cout << "strategy used = "<< clust_seq->strategy_string()<< endl;
652  if (iev == 0) cout << "Jet Definition: " << jet_def.description() << " (" << fj::fastjet_version_string() << ")" << endl;
653  if (do_areas && iev == 0) cout << "Area definition: " << area_def.description() << endl;
654 
655  // now provide some nice output...
656  if (inclkt >= 0.0) {
657  vector<fj::PseudoJet> jets_local = sorted_by_pt(clust_seq->inclusive_jets(inclkt));
658  print_jets(jets_local, show_constituents);
659 
660  }
661 
662  if (excln > 0) {
663  cout << "Printing "<<excln<<" exclusive jets\n";
664  print_jets(clust_seq->exclusive_jets(excln), show_constituents);
665  }
666 
667  if (excld > 0.0) {
668  cout << "Printing exclusive jets for d = "<<excld<<"\n";
669  print_jets(clust_seq->exclusive_jets(excld), show_constituents);
670  }
671 
672  if (excly > 0.0) {
673  cout << "Printing exclusive jets for ycut = "<<excly<<"\n";
674  print_jets(clust_seq->exclusive_jets_ycut(excly), show_constituents);
675  }
676 
677  if (get_all_dij) {
678  for (int i = nparticles-1; i >= 0; i--) {
679  printf("d for n = %4d -> %4d is %14.5e\n", i+1, i, clust_seq->exclusive_dmerge(i));
680  }
681  }
682  if (get_all_yij) {
683  for (int i = nparticles-1; i >= 0; i--) {
684  printf("y for n = %4d -> %4d is %14.5e\n", i+1, i, clust_seq->exclusive_ymerge(i));
685  }
686  }
687 
688  // have the option of printing out the subjets (at scale dcut) of
689  // each inclusive jet
690  if (subdcut >= 0.0) {
691  print_jets_and_sub(clust_seq->inclusive_jets(), subdcut);
692  }
693 
694  // useful for testing that recombination sequences are unique
695  if (unique_write) {
696  vector<int> unique_history = clust_seq->unique_history_order();
697  // construct the inverse of the above mapping
698  vector<int> inv_unique_history(clust_seq->history().size());
699  for (unsigned int i = 0; i < unique_history.size(); i++) {
700  inv_unique_history[unique_history[i]] = i;}
701 
702  for (unsigned int i = 0; i < unique_history.size(); i++) {
703  fj::ClusterSequence::history_element el =
704  clust_seq->history()[unique_history[i]];
705  int uhp1 = el.parent1>=0 ? inv_unique_history[el.parent1] : el.parent1;
706  int uhp2 = el.parent2>=0 ? inv_unique_history[el.parent2] : el.parent2;
707  printf("%7d u %15.8e %7d u %7d u\n",i,el.dij,uhp1, uhp2);
708  }
709  }
710 
711 
712 #ifdef FASTJET_ENABLE_PLUGIN_SISCONE
713  // provide some complementary information for SISCone
714  if (show_cones) {
715  const fj::SISConeExtras * extras =
716  dynamic_cast<const fj::SISConeExtras *>(clust_seq->extras());
717  if (extras == 0)
718  throw fastjet::Error("extras object for SISCone was null (this can happen with certain area types)");
719  cout << "most ambiguous split (difference in squared dist) = "
720  << extras->most_ambiguous_split() << endl;
721  vector<fastjet::PseudoJet> stable_cones(extras->stable_cones());
722  stable_cones = sorted_by_rapidity(stable_cones);
723  for (unsigned int i = 0; i < stable_cones.size(); i++) {
724  //if (stable_cones[i].phi() < 5.0 && stable_cones[i].phi() > 4.0) {
725  printf("%5u %15.8f %15.8f %15.8f\n",
726  i,stable_cones[i].rap(),stable_cones[i].phi(),
727  stable_cones[i].perp() );
728  //}
729  }
730 
731  // also show passes for jets
732  vector<fj::PseudoJet> sisjets = clust_seq->inclusive_jets();
733  printf("\n%15s %15s %15s %12s %8s %8s\n","rap","phi","pt","user-index","pass","nconst");
734  for (unsigned i = 0; i < sisjets.size(); i++) {
735  printf("%15.8f %15.8f %15.8f %12d %8d %8u\n",
736  sisjets[i].rap(), sisjets[i].phi(), sisjets[i].perp(),
737  sisjets[i].user_index(), extras->pass(sisjets[i]),
738  (unsigned int) clust_seq->constituents(sisjets[i]).size()
739  );
740 
741  }
742  }
743 #endif // FASTJET_ENABLE_PLUGIN_SISCONE
744 
745  if (do_bkgd) {
746  double rho, sigma, mean_area, empty_area, n_empty_jets;
748  dynamic_cast<fj::ClusterSequenceAreaBase *>(clust_seq.get());
749  if (do_bkgd_csab) {
750  csab->get_median_rho_and_sigma(bkgd_range, true, rho, sigma, mean_area);
751  empty_area = csab->empty_area(bkgd_range);
752  n_empty_jets = csab->n_empty_jets(bkgd_range);
753  } else if (do_bkgd_jetmedian) {
754  fj::JetMedianBackgroundEstimator bge(bkgd_range);
755  bge.set_provide_fj2_sigma(do_bkgd_fj2);
756  bge.set_cluster_sequence(*csab);
757  rho = bge.rho();
758  sigma = bge.sigma();
759  mean_area = bge.mean_area();
760  empty_area = bge.empty_area();
761  n_empty_jets = bge.n_empty_jets();
762  } else {
763  assert(do_bkgd_gridmedian);
764  double rapmin, rapmax;
765  bkgd_range.get_rapidity_extent(rapmin, rapmax);
766  fj::GridMedianBackgroundEstimator bge(rapmax, 2*ktR);
767  bge.set_particles(particles);
768  rho = bge.rho();
769  sigma = bge.sigma();
770  mean_area = bge.mean_area();
771  empty_area = 0;
772  n_empty_jets = 0;
773  }
774  cout << " rho = " << rho
775  << ", sigma = " << sigma
776  << ", mean_area = " << mean_area
777  << ", empty_area = " << empty_area
778  << ", n_empty_jets = " << n_empty_jets
779  << endl;
780  }
781  } // try
782  catch (fastjet::Error fjerr) {
783  cout << "Caught fastjet error, exiting gracefully" << endl;
784  exit(0);
785  }
786 
787  } // irepeat
788  } // iev
789  // if we've instantiated a plugin, delete it
790  if (jet_def.strategy()==fj::plugin_strategy){
791  delete jet_def.plugin();
792  }
793  // close any file that we've opened
794  if (istr != &cin) delete istr;
795  } // jet_defs
796 
797 }
798 
799 
800 
801 
802 //------ HELPER ROUTINES -----------------------------------------------
803 /// print a single jet
804 void print_jet (const fj::PseudoJet & jet) {
805  unsigned int n_constituents = jet.constituents().size();
806  printf("%15.8f %15.8f %15.8f %8u\n",
807  jet.rap(), jet.phi(), jet.perp(), n_constituents);
808 }
809 
810 
811 //----------------------------------------------------------------------
812 void print_jets(const vector<fj::PseudoJet> & jets_in, bool show_constituents) {
813  vector<fj::PseudoJet> jets;
814  if (ee_print) {
815  jets = sorted_by_E(jets_in);
816  for (unsigned int j = 0; j < jets.size(); j++) {
817  printf("%5u %15.8f %15.8f %15.8f %15.8f\n",
818  j,jets[j].px(),jets[j].py(),jets[j].pz(),jets[j].E());
819  if (show_constituents) {
820  vector<fj::PseudoJet> const_jets = jets[j].constituents();
821  for (unsigned int k = 0; k < const_jets.size(); k++) {
822  printf(" jet%03u %15.8f %15.8f %15.8f %15.8f\n",j,const_jets[k].px(),
823  const_jets[k].py(),const_jets[k].pz(),const_jets[k].E());
824  }
825  cout << "\n\n";
826  }
827 
828  }
829  } else {
830  jets = sorted_by_pt(jets_in);
831  for (unsigned int j = 0; j < jets.size(); j++) {
832  printf("%5u %15.8f %15.8f %15.8f",
833  j,jets[j].rap(),jets[j].phi(),jets[j].perp());
834  // also print out the scalar area and the perp component of the
835  // 4-vector (just enough to check a reasonable 4-vector?)
836  if (do_areas) printf(" %15.8f %15.8f", jets[j].area(),
837  jets[j].area_4vector().perp());
838  cout << "\n";
839 
840  if (show_constituents) {
841  vector<fj::PseudoJet> const_jets = jets[j].constituents();
842  for (unsigned int k = 0; k < const_jets.size(); k++) {
843  printf(" jet%03u %15.8f %15.8f %15.8f %5d\n",j,const_jets[k].rap(),
844  const_jets[k].phi(),sqrt(const_jets[k].kt2()), const_jets[k].cluster_hist_index());
845  }
846  cout << "\n\n";
847  }
848  }
849  }
850 
851  if (rootfile != "") {
852  ofstream ostr(rootfile.c_str());
853  ostr << "# " << cmdline_p->command_line() << endl;
854  ostr << "# output for root" << endl;
855  assert(jets.size() > 0);
856  jets[0].validated_cs()->print_jets_for_root(jets,ostr);
857  }
858 
859 }
860 
861 
862 //----- SUBJETS --------------------------------------------------------
863 /// a function that pretty prints a list of jets and the subjets for each
864 /// one
865 void print_jets_and_sub (const vector<fj::PseudoJet> & jets, double dcut) {
866 
867  // sort jets into increasing pt
868  vector<fj::PseudoJet> sorted_jets = sorted_by_pt(jets);
869 
870  // label the columns
871  printf("Printing jets and their subjets with subdcut = %10.5f\n",dcut);
872  printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity",
873  "phi", "pt", "n constituents");
874 
875  // have various kinds of subjet finding, to test consistency among them
876  enum SubType {internal, newclust_dcut, newclust_R};
877  SubType subtype = internal;
878  //SubType subtype = newclust_dcut;
879  //SubType subtype = newclust_R;
880 
881  // print out the details for each jet
882  //for (unsigned int i = 0; i < sorted_jets.size(); i++) {
883  for (vector<fj::PseudoJet>::const_iterator jet = sorted_jets.begin();
884  jet != sorted_jets.end(); jet++) {
885  const fj::JetDefinition & jet_def = jet->validated_cs()->jet_def();
886 
887  // if jet pt^2 < dcut with kt alg, then some methods of
888  // getting subjets will return nothing -- so skip the jet
889  if (jet_def.jet_algorithm() == fj::kt_algorithm
890  && jet->perp2() < dcut) continue;
891 
892 
893  printf("%5u ",(unsigned int) (jet - sorted_jets.begin()));
894  print_jet(*jet);
895  vector<fj::PseudoJet> subjets;
896  fj::ClusterSequence * cspoint;
897  if (subtype == internal) {
898  cspoint = 0;
899  subjets = jet->exclusive_subjets(dcut);
900  double ddnp1 = jet->exclusive_subdmerge_max(subjets.size());
901  double ddn = jet->exclusive_subdmerge_max(subjets.size()-1);
902  cout << " for " << ddnp1 << " < d < " << ddn << " one has " << endl;
903  } else if (subtype == newclust_dcut) {
904  cspoint = new fj::ClusterSequence(jet->constituents(), jet_def);
905  subjets = cspoint->exclusive_jets(dcut);
906  } else if (subtype == newclust_R) {
907  assert(jet_def.jet_algorithm() == fj::cambridge_algorithm);
908  fj::JetDefinition subjd(jet_def.jet_algorithm(), jet_def.R()*sqrt(dcut));
909  cspoint = new fj::ClusterSequence(jet->constituents(), subjd);
910  subjets = cspoint->inclusive_jets();
911  } else {
912  cerr << "unrecognized subtype for subjet finding" << endl;
913  exit(-1);
914  }
915 
916  subjets = sorted_by_pt(subjets);
917  for (unsigned int j = 0; j < subjets.size(); j++) {
918  printf(" -sub-%02u ",j);
919  print_jet(subjets[j]);
920  }
921 
922  if (cspoint != 0) delete cspoint;
923 
924  //fj::ClusterSequence subseq(clust_seq->constituents(sorted_jets[i]),
925  // fj::JetDefinition(fj::cambridge_algorithm, 0.4));
926  //vector<fj::PseudoJet> subjets = sorted_by_pt(subseq.inclusive_jets());
927  //for (unsigned int j = 0; j < subjets.size(); j++) {
928  // printf(" -sub-%02u ",j);
929  // print_jet(subseq, subjets[j]);
930  //}
931  }
932 
933 }
934 
fastjet::SISConePlugin
Definition: SISConePlugin.hh:71
fastjet::ClusterSequenceAreaBase::n_empty_jets
virtual double n_empty_jets(const Selector &selector) const
return something similar to the number of pure ghost jets in the given selector's range in an active ...
Definition: ClusterSequenceAreaBase.hh:133
fastjet::SelectorAbsRapMax
Selector SelectorAbsRapMax(double absrapmax)
select objects with |rap| <= absrapmax
Definition: Selector.cc:828
fastjet::EECambridgePlugin
Definition: EECambridgePlugin.hh:57
fastjet
Definition: AreaDefinition.cc:35
fastjet::SelectorNHardest
Selector SelectorNHardest(unsigned int n)
select the n hardest objects
Definition: Selector.cc:1022
fastjet::PseudoJet::exclusive_subjets
std::vector< PseudoJet > exclusive_subjets(const double &dcut) const
return a vector of all subjets of the current jet (in the sense of the exclusive algorithm) that woul...
Definition: PseudoJet.cc:589
fastjet::E_scheme
@ E_scheme
summing the 4-momenta
Definition: JetDefinition.hh:136
fastjet::sorted_by_rapidity
vector< PseudoJet > sorted_by_rapidity(const vector< PseudoJet > &jets)
return a vector of jets sorted into increasing rapidity
Definition: PseudoJet.cc:793
fastjet::JetDefinition
Definition: JetDefinition.hh:173
fastjet::PseudoJet::exclusive_subdmerge_max
double exclusive_subdmerge_max(int nsub) const
return the maximum dij that occurred in the whole event at the stage that the nsub+1 -> nsub merge of...
Definition: PseudoJet.cc:648
fastjet::SISConeSphericalPlugin
Definition: SISConeSphericalPlugin.hh:94
fastjet::Selector::get_rapidity_extent
void get_rapidity_extent(double &rapmin, double &rapmax) const
returns the rapidity range for which it may return "true"
Definition: Selector.hh:213
fastjet::VoronoiAreaSpec
Definition: AreaDefinition.hh:48
fastjet::PseudoJet::constituents
virtual std::vector< PseudoJet > constituents() const
retrieve the constituents.
Definition: PseudoJet.cc:566
fastjet::ClusterSequence
Definition: ClusterSequence.hh:59
fastjet::JetMedianBackgroundEstimator
Definition: JetMedianBackgroundEstimator.hh:78
fastjet::ClusterSequence::jet_def
const JetDefinition & jet_def() const
return a reference to the jet definition
Definition: ClusterSequence.hh:291
fastjet::D0RunIIConePlugin
Definition: D0RunIIConePlugin.hh:58
fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma
virtual void get_median_rho_and_sigma(const Selector &selector, bool use_area_4vector, double &median, double &sigma, double &mean_area) const
using jets withing the selector range (and with 4-vector areas if use_area_4vector),...
Definition: ClusterSequenceAreaBase.cc:162
fastjet::sorted_by_E
vector< PseudoJet > sorted_by_E(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing energy
Definition: PseudoJet.cc:801
print_jets
void print_jets(const vector< fastjet::PseudoJet > &)
a function that pretty prints a list of jets
Definition: fastjet_example.cc:126
fastjet::PseudoJet::perp
double perp() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:141
fastjet::JetDefinition::description
std::string description() const
return a textual description of the current jet definition
Definition: JetDefinition.cc:108
fastjet::PseudoJet::validated_cs
const ClusterSequence * validated_cs() const
shorthand for validated_cluster_sequence()
Definition: PseudoJet.cc:443
fastjet::TrackJetPlugin
Definition: TrackJetPlugin.hh:46
fastjet::AreaDefinition
Definition: AreaDefinition.hh:80
fastjet::ClusterSequenceArea
Definition: ClusterSequenceArea.hh:49
fastjet::ClusterSequenceAreaBase::empty_area
virtual double empty_area(const Selector &selector) const
return the total area, corresponding to the given Selector, that is free of jets, in general based on...
Definition: ClusterSequenceAreaBase.cc:55
fastjet::fastjet_version_string
string fastjet_version_string()
return a string containing information about the release
Definition: ClusterSequence.cc:366
fastjet::JetDefinition::jet_algorithm
JetAlgorithm jet_algorithm() const
return information about the definition...
Definition: JetDefinition.hh:317
fastjet::SISConeBaseExtras::pass
int pass(const PseudoJet &jet) const
return the # of the pass at which a given jet was found; will return -1 if the pass is invalid
Definition: SISConeBasePlugin.hh:161
fastjet::PxConePlugin
Definition: PxConePlugin.hh:66
main
int main()
an example program showing how to use fastjet
Definition: 01-basic.cc:50
fastjet::ClusterSequence::print_banner
static void print_banner()
This is the function that is automatically called during clustering to print the FastJet banner.
Definition: ClusterSequence.cc:373
fastjet::GhostedAreaSpec
Definition: GhostedAreaSpec.hh:62
fastjet::ClusterSequence::inclusive_jets
std::vector< PseudoJet > inclusive_jets(const double &ptmin=0.0) const
return a vector of all jets (in the sense of the inclusive algorithm) with pt >= ptmin.
Definition: ClusterSequence.cc:664
fastjet::PseudoJet
Definition: PseudoJet.hh:65
fastjet::ClusterSequenceAreaBase
Definition: ClusterSequenceAreaBase.hh:45
fastjet::PseudoJet::rap
double rap() const
returns the rapidity or some large value when the rapidity is infinite
Definition: PseudoJet.hh:121
fastjet::ClusterSequence::exclusive_jets
std::vector< PseudoJet > exclusive_jets(const double &dcut) const
return a vector of all jets (in the sense of the exclusive algorithm) that would be obtained when run...
Definition: ClusterSequence.cc:735
fastjet::Strategy
Strategy
the various options for the algorithmic strategy to adopt in clustering events with kt and cambridge ...
Definition: JetDefinition.hh:48
fastjet::Selector
Definition: Selector.hh:141
fastjet::PseudoJet::phi
double phi() const
returns phi (in the range 0..2pi)
Definition: PseudoJet.hh:106
fastjet::JadePlugin
Definition: JadePlugin.hh:75
fastjet::SISConeBaseExtras::stable_cones
const std::vector< PseudoJet > & stable_cones() const
returns a reference to the vector of stable cones (aka protocones)
Definition: SISConeBasePlugin.hh:154
fastjet::RecombinationScheme
RecombinationScheme
the various recombination schemes
Definition: JetDefinition.hh:134
fastjet::PseudoJet::perp2
double perp2() const
returns the squared transverse momentum
Definition: PseudoJet.hh:139
fastjet::CMSIterativeConePlugin
Definition: CMSIterativeConePlugin.hh:46
fastjet::JetDefinition::plugin
const Plugin * plugin() const
return a pointer to the plugin
Definition: JetDefinition.hh:310
fastjet::sorted_by_pt
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
Definition: PseudoJet.cc:785
fastjet::GridJetPlugin
Definition: GridJetPlugin.hh:51
fastjet::Error
Definition: Error.hh:41
fastjet::SISConeBaseExtras::most_ambiguous_split
double most_ambiguous_split() const
return the smallest difference in squared distance encountered during splitting between a particle an...
Definition: SISConeBasePlugin.hh:175
fastjet::SISConeExtras
Definition: SISConePlugin.hh:198