FastJet 3.0.4
JetDefinition.hh
00001 //STARTHEADER
00002 // $Id: JetDefinition.hh 2687 2011-11-14 11:17:51Z soyez $
00003 //
00004 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
00005 //
00006 //----------------------------------------------------------------------
00007 // This file is part of FastJet.
00008 //
00009 //  FastJet is free software; you can redistribute it and/or modify
00010 //  it under the terms of the GNU General Public License as published by
00011 //  the Free Software Foundation; either version 2 of the License, or
00012 //  (at your option) any later version.
00013 //
00014 //  The algorithms that underlie FastJet have required considerable
00015 //  development and are described in hep-ph/0512210. If you use
00016 //  FastJet as part of work towards a scientific publication, please
00017 //  include a citation to the FastJet paper.
00018 //
00019 //  FastJet is distributed in the hope that it will be useful,
00020 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00021 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00022 //  GNU General Public License for more details.
00023 //
00024 //  You should have received a copy of the GNU General Public License
00025 //  along with FastJet. If not, see <http://www.gnu.org/licenses/>.
00026 //----------------------------------------------------------------------
00027 //ENDHEADER
00028 
00029 #ifndef __FASTJET_JETDEFINITION_HH__
00030 #define __FASTJET_JETDEFINITION_HH__
00031 
00032 #include<cassert>
00033 #include "fastjet/internal/numconsts.hh"
00034 #include "fastjet/PseudoJet.hh"
00035 #include<string>
00036 #include<memory>
00037 
00038 FASTJET_BEGIN_NAMESPACE      // defined in fastjet/internal/base.hh
00039 
00040 /// return a string containing information about the release
00041 //  NB: (implemented in ClusterSequence.cc but defined here because
00042 //  this is a visible location)
00043 std::string fastjet_version_string();
00044 
00045 //======================================================================
00046 /// the various options for the algorithmic strategy to adopt in
00047 /// clustering events with kt and cambridge style algorithms.
00048 enum Strategy {
00049   /// fastest form about 500..10^4
00050   N2MinHeapTiled   = -4, 
00051   /// fastest from about 50..500
00052   N2Tiled     = -3, 
00053   /// legacy
00054   N2PoorTiled = -2, 
00055   /// fastest below 50
00056   N2Plain     = -1, 
00057   /// worse even than the usual N^3 algorithms
00058   N3Dumb      =  0, 
00059   /// automatic selection of the best (based on N)
00060   Best        =  1, 
00061   /// best of the NlnN variants -- best overall for N>10^4.
00062   /// (Does not work for R>=2pi)
00063   NlnN        =  2, 
00064   /// legacy N ln N using 3pi coverage of cylinder.
00065   /// (Does not work for R>=2pi)
00066   NlnN3pi     =  3, 
00067   /// legacy N ln N using 4pi coverage of cylinder
00068   NlnN4pi     =  4,
00069   /// Chan's closest pair method (in a variant with 4pi coverage),
00070   /// for use exclusively with the Cambridge algorithm.
00071   /// (Does not work for R>=2pi)
00072   NlnNCam4pi   = 14,
00073   /// Chan's closest pair method (in a variant with 2pi+2R coverage),
00074   /// for use exclusively with the Cambridge algorithm.
00075   /// (Does not work for R>=2pi)
00076   NlnNCam2pi2R = 13,
00077   /// Chan's closest pair method (in a variant with 2pi+minimal extra
00078   /// variant), for use exclusively with the Cambridge algorithm. 
00079   /// (Does not work for R>=2pi)
00080   NlnNCam      = 12, // 2piMultD
00081   /// the plugin has been used...
00082   plugin_strategy = 999
00083 };
00084 
00085 
00086 //======================================================================
00087 /// \enum JetAlgorithm
00088 /// the various families of jet-clustering algorithm
00089 enum JetAlgorithm {
00090   /// the longitudinally invariant kt algorithm
00091   kt_algorithm=0,
00092   /// the longitudinally invariant variant of the cambridge algorithm
00093   /// (aka Aachen algoithm).
00094   cambridge_algorithm=1,
00095   /// like the k_t but with distance measures 
00096   ///       dij = min(1/kti^2,1/ktj^2) Delta R_{ij}^2 / R^2
00097   ///       diB = 1/kti^2
00098   antikt_algorithm=2, 
00099   /// like the k_t but with distance measures 
00100   ///       dij = min(kti^{2p},ktj^{2p}) Delta R_{ij}^2 / R^2
00101   ///       diB = 1/kti^{2p}
00102   /// where p = extra_param()
00103   genkt_algorithm=3, 
00104   /// a version of cambridge with a special distance measure for particles
00105   /// whose pt is < extra_param()
00106   cambridge_for_passive_algorithm=11,
00107   /// a version of genkt with a special distance measure for particles
00108   /// whose pt is < extra_param() [relevant for passive areas when p<=0]
00109   genkt_for_passive_algorithm=13, 
00110   //.................................................................
00111   /// the e+e- kt algorithm
00112   ee_kt_algorithm=50,
00113   /// the e+e- genkt algorithm  (R > 2 and p=1 gives ee_kt)
00114   ee_genkt_algorithm=53,
00115   //.................................................................
00116   /// any plugin algorithm supplied by the user
00117   plugin_algorithm = 99,
00118   //.................................................................
00119   /// the value for the jet algorithm in a JetDefinition for which
00120   /// no algorithm has yet been defined
00121   undefined_jet_algorithm = 999
00122 };
00123 
00124 /// make standard Les Houches nomenclature JetAlgorithm (algorithm is general
00125 /// recipe without the parameters) backward-compatible with old JetFinder
00126 typedef JetAlgorithm JetFinder;
00127 
00128 /// provide other possible names for the Cambridge/Aachen algorithm
00129 const JetAlgorithm aachen_algorithm = cambridge_algorithm;
00130 const JetAlgorithm cambridge_aachen_algorithm = cambridge_algorithm;
00131 
00132 //======================================================================
00133 /// the various recombination schemes
00134 enum RecombinationScheme {
00135   /// summing the 4-momenta
00136   E_scheme=0,
00137   /// pt weighted recombination of y,phi (and summing of pt's)
00138   /// with preprocessing to make things massless by rescaling E=|\vec p|
00139   pt_scheme=1,
00140   /// pt^2 weighted recombination of y,phi (and summing of pt's)
00141   /// with preprocessing to make things massless by rescaling E=|\vec p|
00142   pt2_scheme=2,
00143   /// pt weighted recombination of y,phi (and summing of pt's)
00144   /// with preprocessing to make things massless by rescaling |\vec p|->=E
00145   Et_scheme=3,
00146   /// pt^2 weighted recombination of y,phi (and summing of pt's)
00147   /// with preprocessing to make things massless by rescaling |\vec p|->=E
00148   Et2_scheme=4,
00149   /// pt weighted recombination of y,phi (and summing of pt's), with 
00150   /// no preprocessing
00151   BIpt_scheme=5,
00152   /// pt^2 weighted recombination of y,phi (and summing of pt's)
00153   /// no preprocessing
00154   BIpt2_scheme=6,
00155   /// for the user's external scheme
00156   external_scheme = 99
00157 };
00158 
00159 
00160 
00161 // forward declaration, needed in order to specify interface for the
00162 // plugin.
00163 class ClusterSequence;
00164 
00165 
00166 
00167 
00168 //======================================================================
00169 /// @ingroup basic_classes
00170 /// \class JetDefinition
00171 /// class that is intended to hold a full definition of the jet
00172 /// clusterer
00173 class JetDefinition {
00174   
00175 public:
00176 
00177   /// forward declaration of a class that allows the user to introduce
00178   /// their own plugin 
00179   class Plugin;
00180 
00181   // forward declaration of a class that will provide the
00182   // recombination scheme facilities and/or allow a user to
00183   // extend these facilities
00184   class Recombiner;
00185 
00186 
00187   /// constructor with alternative ordering or arguments -- note that
00188   /// we have not provided a default jet finder, to avoid ambiguous
00189   /// JetDefinition() constructor.
00190   JetDefinition(JetAlgorithm jet_algorithm_in, 
00191                 double R_in, 
00192                 RecombinationScheme recomb_scheme_in = E_scheme,
00193                 Strategy strategy_in = Best) {
00194     *this = JetDefinition(jet_algorithm_in, R_in, strategy_in, recomb_scheme_in, 1);
00195   }
00196 
00197   /// constructor for algorithms that have no free parameters
00198   /// (e.g. ee_kt_algorithm)
00199   JetDefinition(JetAlgorithm jet_algorithm_in, 
00200                 RecombinationScheme recomb_scheme_in = E_scheme,
00201                 Strategy strategy_in = Best) {
00202     double dummyR = 0.0;
00203     *this = JetDefinition(jet_algorithm_in, dummyR, strategy_in, recomb_scheme_in, 0);
00204   }
00205 
00206   /// constructor for algorithms that require R + one extra parameter to be set 
00207   /// (the gen-kt series for example)
00208   JetDefinition(JetAlgorithm jet_algorithm_in, 
00209                 double R_in, 
00210                 double xtra_param_in,
00211                 RecombinationScheme recomb_scheme_in = E_scheme,
00212                 Strategy strategy_in = Best) {
00213     *this = JetDefinition(jet_algorithm_in, R_in, strategy_in, recomb_scheme_in, 2);
00214     set_extra_param(xtra_param_in);
00215   }
00216 
00217 
00218   /// constructor in a form that allows the user to provide a pointer
00219   /// to an external recombiner class (which must remain valid for the
00220   /// life of the JetDefinition object).
00221   JetDefinition(JetAlgorithm jet_algorithm_in, 
00222                 double R_in, 
00223                 const Recombiner * recombiner_in,
00224                 Strategy strategy_in = Best) {
00225     *this = JetDefinition(jet_algorithm_in, R_in, external_scheme, strategy_in);
00226     _recombiner = recombiner_in;
00227   }
00228 
00229 
00230   /// constructor for case with 0 parameters (ee_kt_algorithm) and
00231   /// and external recombiner
00232   JetDefinition(JetAlgorithm jet_algorithm_in, 
00233                 const Recombiner * recombiner_in,
00234                 Strategy strategy_in = Best) {
00235     *this = JetDefinition(jet_algorithm_in, external_scheme, strategy_in);
00236     _recombiner = recombiner_in;
00237   }
00238 
00239   /// constructor allowing the extra parameter to be set and a pointer to
00240   /// a recombiner
00241   JetDefinition(JetAlgorithm jet_algorithm_in, 
00242                 double R_in, 
00243                 double xtra_param_in,
00244                 const Recombiner * recombiner_in,
00245                 Strategy strategy_in = Best) {
00246     *this = JetDefinition(jet_algorithm_in, R_in, external_scheme, strategy_in);
00247     _recombiner = recombiner_in;
00248     set_extra_param(xtra_param_in);
00249   }
00250 
00251   /// a default constructor which creates a jet definition that is in
00252   /// a well-defined internal state, but not actually usable for jet
00253   /// clustering.
00254   JetDefinition()  {
00255     *this = JetDefinition(undefined_jet_algorithm, 1.0);
00256   }
00257   
00258 
00259   // /// a default constructor
00260   // JetDefinition() {
00261   //   *this = JetDefinition(kt_algorithm, 1.0);
00262   // }
00263 
00264   /// constructor based on a pointer to a user's plugin; the object
00265   /// pointed to must remain valid for the whole duration of existence
00266   /// of the JetDefinition and any related ClusterSequences
00267   JetDefinition(const Plugin * plugin_in) {
00268     _plugin = plugin_in;
00269     _strategy = plugin_strategy;
00270     _Rparam = _plugin->R();
00271     _jet_algorithm = plugin_algorithm;
00272     set_recombination_scheme(E_scheme);
00273   }
00274 
00275 
00276   /// constructor to fully specify a jet-definition (together with
00277   /// information about how algorithically to run it).
00278   ///
00279   /// the ordering of arguments here is old and deprecated (except
00280   /// as the common constructor for internal use)
00281   JetDefinition(JetAlgorithm jet_algorithm_in, 
00282                 double R_in, 
00283                 Strategy strategy_in,
00284                 RecombinationScheme recomb_scheme_in = E_scheme,
00285                 int nparameters_in = 1);
00286   
00287   /// R values larger than max_allowable_R are not allowed.
00288   ///
00289   /// We use a value of 1000, substantially smaller than
00290   /// numeric_limits<double>::max(), to leave room for the convention
00291   /// within PseudoJet of setting unphysical (infinite) rapidities to
00292   /// +-(MaxRap + abs(pz())), where MaxRap is 10^5.
00293   static const double max_allowable_R; //= 1000.0;
00294 
00295   /// set the recombination scheme to the one provided
00296   void set_recombination_scheme(RecombinationScheme);
00297 
00298   /// set the recombiner class to the one provided
00299   void set_recombiner(const Recombiner * recomb) {
00300     if (_recombiner_shared()) _recombiner_shared.reset(recomb);
00301     _recombiner = recomb;
00302     _default_recombiner = DefaultRecombiner(external_scheme);
00303   }
00304 
00305   /// calling this tells the JetDefinition to handle the deletion of
00306   /// the recombiner when it is no longer used
00307   void delete_recombiner_when_unused();
00308 
00309   /// return a pointer to the plugin 
00310   const Plugin * plugin() const {return _plugin;};
00311 
00312   /// allows to let the JetDefinition handle the deletion of the
00313   /// plugin when it is no longer used
00314   void delete_plugin_when_unused();
00315 
00316   /// return information about the definition...
00317   JetAlgorithm jet_algorithm  () const {return _jet_algorithm  ;}
00318   /// same as above for backward compatibility
00319   JetAlgorithm jet_finder     () const {return _jet_algorithm  ;}
00320   double    R           () const {return _Rparam      ;}
00321   // a general purpose extra parameter, whose meaning depends on
00322   // the algorithm, and may often be unused.
00323   double    extra_param () const {return _extra_param ;}
00324   Strategy  strategy    () const {return _strategy    ;}
00325   RecombinationScheme recombination_scheme() const {
00326     return _default_recombiner.scheme();}
00327 
00328   /// (re)set the jet finder
00329   void set_jet_algorithm(JetAlgorithm njf) {_jet_algorithm = njf;}
00330   /// same as above for backward compatibility
00331   void set_jet_finder(JetAlgorithm njf)    {_jet_algorithm = njf;}
00332   /// (re)set the general purpose extra parameter
00333   void set_extra_param(double xtra_param) {_extra_param = xtra_param;}
00334 
00335   /// return a pointer to the currently defined recombiner. 
00336   ///
00337   /// Warning: the pointer may be to an internal recombiner (for
00338   /// default recombination schemes), in which case if the
00339   /// JetDefinition becomes invalid (e.g. is deleted), the pointer
00340   /// will then point to an object that no longer exists.
00341   /// 
00342   /// Note also that if you copy a JetDefinition with a default
00343   /// recombination scheme, then the two copies will have distinct
00344   /// recombiners, and return different recombiner() pointers.
00345   const Recombiner * recombiner() const {
00346     return _recombiner == 0 ? & _default_recombiner : _recombiner;}
00347 
00348   /// returns true if the current jet definitions shares the same
00349   /// recombiner as teh one passed as an argument
00350   bool has_same_recombiner(const JetDefinition &other_jd) const;
00351 
00352   /// return a textual description of the current jet definition 
00353   std::string description() const;
00354 
00355 
00356 public:
00357   //======================================================================
00358   /// @ingroup advanced_usage
00359   /// \class Recombiner
00360   /// An abstract base class that will provide the recombination scheme
00361   /// facilities and/or allow a user to extend these facilities
00362   class Recombiner {
00363   public:
00364     /// return a textual description of the recombination scheme
00365     /// implemented here
00366     virtual std::string description() const = 0;
00367     
00368     /// recombine pa and pb and put result into pab
00369     virtual void recombine(const PseudoJet & pa, const PseudoJet & pb, 
00370                            PseudoJet & pab) const = 0;
00371 
00372     /// routine called to preprocess each input jet (to make all input
00373     /// jets compatible with the scheme requirements (e.g. massless).
00374     virtual void preprocess(PseudoJet & ) const {};
00375     
00376     /// a destructor to be replaced if necessary in derived classes...
00377     virtual ~Recombiner() {};
00378 
00379     /// pa += pb in the given recombination scheme. Not virtual -- the
00380     /// user should have no reason to want to redefine this!
00381     inline void plus_equal(PseudoJet & pa, const PseudoJet & pb) const {
00382       // put result in a temporary location in case the recombiner
00383       // does something funny (ours doesn't, but who knows about the
00384       // user's)
00385       PseudoJet pres; 
00386       recombine(pa,pb,pres);
00387       pa = pres;
00388     }
00389 
00390   };
00391   
00392   
00393   //======================================================================
00394   /// @ingroup advanced_usage
00395   /// \class DefaultRecombiner
00396   /// A class that will provide the recombination scheme facilities and/or
00397   /// allow a user to extend these facilities
00398   ///
00399   /// This class is derived from the (abstract) class Recombiner. It
00400   /// simply "sums" PseudoJets using a specified recombination scheme
00401   /// (E-scheme by default)
00402   class DefaultRecombiner : public Recombiner {
00403   public:
00404     DefaultRecombiner(RecombinationScheme recomb_scheme = E_scheme) : 
00405       _recomb_scheme(recomb_scheme) {}
00406     
00407     virtual std::string description() const;
00408     
00409     /// recombine pa and pb and put result into pab
00410     virtual void recombine(const PseudoJet & pa, const PseudoJet & pb, 
00411                            PseudoJet & pab) const;
00412 
00413     virtual void preprocess(PseudoJet & p) const;
00414 
00415     /// return the index of the recombination scheme
00416     RecombinationScheme scheme() const {return _recomb_scheme;}
00417     
00418   private:
00419     RecombinationScheme _recomb_scheme;
00420   };
00421 
00422 
00423   //======================================================================
00424   /// @ingroup advanced_usage
00425   /// \class Plugin
00426   /// a class that allows a user to introduce their own "plugin" jet
00427   /// finder
00428   ///
00429   /// Note that all the plugins provided with FastJet are derived from
00430   /// this class
00431   class Plugin{
00432   public:
00433     /// return a textual description of the jet-definition implemented
00434     /// in this plugin
00435     virtual std::string description() const = 0;
00436     
00437     /// given a ClusterSequence that has been filled up with initial
00438     /// particles, the following function should fill up the rest of the
00439     /// ClusterSequence, using the following member functions of
00440     /// ClusterSequence:
00441     ///   - plugin_do_ij_recombination(...)
00442     ///   - plugin_do_iB_recombination(...)
00443     virtual void run_clustering(ClusterSequence &) const = 0;
00444     
00445     virtual double R() const = 0;
00446     
00447     /// return true if there is specific support for the measurement
00448     /// of passive areas, in the sense that areas determined from all
00449     /// particles below the ghost separation scale will be a passive
00450     /// area. [If you don't understand this, ignore it!]
00451     virtual bool supports_ghosted_passive_areas() const {return false;}
00452 
00453     /// set the ghost separation scale for passive area determinations
00454     /// in future runs (strictly speaking that makes the routine
00455     /// a non const, so related internal info must be stored as a mutable)
00456     virtual void set_ghost_separation_scale(double scale) const;
00457     virtual double ghost_separation_scale() const {return 0.0;}
00458 
00459     /// if this returns false then a warning will be given
00460     /// whenever the user requests "exclusive" jets from the
00461     /// cluster sequence
00462     virtual bool exclusive_sequence_meaningful() const {return false;}
00463 
00464     /// a destructor to be replaced if necessary in derived classes...
00465     virtual ~Plugin() {};
00466   };
00467 
00468 private:
00469 
00470 
00471   JetAlgorithm _jet_algorithm;
00472   double    _Rparam;
00473   double    _extra_param ; ///< parameter whose meaning varies according to context
00474   Strategy  _strategy  ;
00475 
00476   const Plugin * _plugin;
00477   SharedPtr<const Plugin> _plugin_shared;
00478 
00479   // when we use our own recombiner it's useful to point to it here
00480   // so that we don't have to worry about deleting it etc...
00481   DefaultRecombiner _default_recombiner;
00482   const Recombiner * _recombiner;
00483   SharedPtr<const Recombiner> _recombiner_shared;
00484 
00485 };
00486 
00487 
00488 //-------------------------------------------------------------------------------
00489 // helper functions to build a jet made of pieces
00490 //
00491 // These functions include an options recombiner used to compute the
00492 // total composite jet momentum
00493 // -------------------------------------------------------------------------------
00494 
00495 /// build a "CompositeJet" from the vector of its pieces
00496 ///
00497 /// In this case, E-scheme recombination is assumed to compute the
00498 /// total momentum
00499 PseudoJet join(const std::vector<PseudoJet> & pieces, const JetDefinition::Recombiner & recombiner);
00500 
00501 /// build a MergedJet from a single PseudoJet
00502 PseudoJet join(const PseudoJet & j1, 
00503                const JetDefinition::Recombiner & recombiner);
00504 
00505 /// build a MergedJet from 2 PseudoJet
00506 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, 
00507                const JetDefinition::Recombiner & recombiner);
00508 
00509 /// build a MergedJet from 3 PseudoJet
00510 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, 
00511                const JetDefinition::Recombiner & recombiner);
00512 
00513 /// build a MergedJet from 4 PseudoJet
00514 PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, const PseudoJet & j4, 
00515                const JetDefinition::Recombiner & recombiner);
00516 
00517 
00518 
00519 
00520 
00521 FASTJET_END_NAMESPACE
00522 
00523 #endif // __FASTJET_JETDEFINITION_HH__
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends