JETAN/fastjet/fastjet/SISConePlugin.hh

   1 #ifndef __SISCONEPLUGIN_HH__
   2 #define __SISCONEPLUGIN_HH__
   3
   4 #include "SISConeBasePlugin.hh"
   5
   6 // forward declaration of the siscone classes we'll need
   7 namespace siscone{
   8   class Csiscone;
   9 }
  10
  11
  12 namespace fastjet {      // defined in fastjet/internal/base.hh
  13
  14 //----------------------------------------------------------------------
  15 //
  16 /// SISConePlugin is a plugin for fastjet (v2.1 upwards) that provides
  17 /// an interface to the seedless infrared safe cone jet finder by
  18 /// Gregory Soyez and Gavin Salam.
  19 ///
  20 /// SISCone uses geometrical techniques to exhaustively consider all
  21 /// possible distinct cones. It then finds out which ones are stable
  22 /// and sends the result to the Tevatron Run-II type split-merge
  23 /// procedure for overlapping cones.
  24 ///
  25 /// Four parameters govern the "physics" of the algorithm:
  26 ///
  27 ///  - the cone_radius (this should be self-explanatory!)
  28 ///
  29 ///  - the overlap_threshold is the parameter which dictates how much
  30 ///    two jets must overlap (pt_overlap/min(pt1,pt2)) if they are to be
  31 ///    merged
  32 ///
  33 ///  - Not all particles are in stable cones in the first round of
  34 ///    searching for stable cones; one can therefore optionally have the
  35 ///    the jet finder carry out additional passes of searching for
  36 ///    stable cones among particles that were in no stable cone in
  37 ///    previous passes --- the maximum number of passes carried out is
  38 ///    n_pass_max. If this is zero then additional passes are carried
  39 ///    out until no new stable cones are found.
  40 ///
  41 ///  - Protojet ptmin: protojets that are below this ptmin
  42 ///    (default = 0) are discarded before each iteration of the
  43 ///    split-merge loop.
  44 ///
  45 /// One parameter governs some internal algorithmic shortcuts:
  46 ///
  47 /// - if "caching" is turned on then the last event clustered by
  48 ///   siscone is stored -- if the current event is identical and the
  49 ///   cone_radius and n_pass_mass are identical, then the only part of
  50 ///   the clustering that needs to be rerun is the split-merge part,
  51 ///   leading to significant speed gains; there is a small (O(N) storage
  52 ///   and speed) penalty for caching, so it should be kept off
  53 ///   (default) if only a single overlap_threshold is used.
  54 ///
  55 /// The final jets can be accessed by requestion the
  56 /// inclusive_jets(...) from the ClusterSequence object. Note that
  57 /// these PseudoJets have their user_index() set to the index of the
  58 /// pass in which they were found (first pass = 0). NB: This does not
  59 /// currently work for jets that consist of a single particle.
  60 ///
  61 /// For further information on the details of the algorithm see the
  62 /// SISCone paper, arXiv:0704.0292 [JHEP 0705:086,2007].
  63 ///
  64 /// For documentation about the implementation, see the
  65 /// siscone/doc/html/index.html file.
  66 //
  67 class SISConePlugin : public SISConeBasePlugin{
  68 public:
  69
  70   /// enum for the different split-merge scale choices;
  71   /// Note that order _must_ be the same as in siscone
  72   enum SplitMergeScale {SM_pt,     ///< transverse momentum (E-scheme), IR unsafe
  73                         SM_Et,     ///< transverse energy (E-scheme), not long. boost invariant
  74                                    ///< original run-II choice [may not be implemented]
  75                         SM_mt,     ///< transverse mass (E-scheme), IR safe except
  76                                    ///< in decays of two identical narrow heavy particles
  77                         SM_pttilde ///< pt-scheme pt = \sum_{i in jet} |p_{ti}|, should
  78                                    ///< be IR safe in all cases
  79   };
  80
  81
  82   /// Main constructor for the SISCone Plugin class.
  83   ///
  84   /// Note: wrt version prior to 2.4 this constructor differs in that a
  85   /// the default value has been removed for overlap_threshold. The
  86   /// former has been removed because the old default of 0.5 was found
  87   /// to be unsuitable in high-noise environments; so the user should
  88   /// now explicitly think about the value for this -- we recommend
  89   /// 0.75.
  90   ///
  91   SISConePlugin (double cone_radius,
  92                  double overlap_threshold,
  93                  int    n_pass_max = 0,
  94                  double protojet_ptmin = 0.0,
  95                  bool   caching = false,
  96                  SplitMergeScale  split_merge_scale = SM_pttilde,
  97                  double split_merge_stopping_scale = 0.0){
  98     _cone_radius           = cone_radius;
  99     _overlap_threshold     = overlap_threshold;
 100     _n_pass_max            = n_pass_max;
 101     _protojet_ptmin        = protojet_ptmin;
 102     _caching               = caching;
 103     _split_merge_scale     = split_merge_scale;
 104     _split_merge_stopping_scale = split_merge_stopping_scale;
 105     _ghost_sep_scale       = 0.0;
 106     _use_pt_weighted_splitting = false;}
 107
 108
 109   /// Backwards compatible constructor for the SISCone Plugin class
 110   SISConePlugin (double cone_radius,
 111                  double overlap_threshold,
 112                  int    n_pass_max,
 113                  double protojet_ptmin,
 114                  bool   caching ,
 115                  bool   split_merge_on_transverse_mass){
 116     _cone_radius           = cone_radius;
 117     _overlap_threshold     = overlap_threshold;
 118     _n_pass_max            = n_pass_max;
 119     _protojet_ptmin        = protojet_ptmin;
 120     _caching               = caching;
 121     _split_merge_stopping_scale = 0.0;
 122     _split_merge_scale     = split_merge_on_transverse_mass ? SM_mt : SM_pttilde;
 123     _ghost_sep_scale       = 0.0;}
 124
 125   /// backwards compatible constructor for the SISCone Plugin class
 126   /// (avoid using this in future).
 127   SISConePlugin (double cone_radius,
 128                  double overlap_threshold,
 129                  int    n_pass_max,
 130                  bool   caching ) {
 131     _cone_radius           = cone_radius;
 132     _overlap_threshold     = overlap_threshold;
 133     _n_pass_max            = n_pass_max;
 134     _protojet_ptmin        = 0.0;
 135     _caching               = caching;
 136     _split_merge_scale     = SM_mt;
 137     _split_merge_stopping_scale = 0.0;
 138     _ghost_sep_scale       = 0.0;
 139     _use_pt_weighted_splitting = false;}
 140
 141   /// minimum pt for a protojet to be considered in the split-merge step
 142   /// of the algorithm
 143   double protojet_ptmin  () const {return _protojet_ptmin  ;}
 144
 145   /// return the scale to be passed to SISCone as the protojet_ptmin
 146   /// -- if we have a ghost separation scale that is above the
 147   /// protojet_ptmin, then the ghost_separation_scale becomes the
 148   /// relevant one to use here
 149   double protojet_or_ghost_ptmin  () const {return std::max(_protojet_ptmin,
 150                                                             _ghost_sep_scale);}
 151
 152   /// indicates scale used in split-merge
 153   SplitMergeScale split_merge_scale() const {return _split_merge_scale;}
 154   /// sets scale used in split-merge
 155   void set_split_merge_scale(SplitMergeScale sms) {_split_merge_scale = sms;}
 156
 157   /// indicates whether the split-merge orders on transverse mass or not.
 158   /// retained for backwards compatibility with 2.1.0b3
 159   bool split_merge_on_transverse_mass() const {return _split_merge_scale == SM_mt ;}
 160   void set_split_merge_on_transverse_mass(bool val) {
 161     _split_merge_scale = val  ? SM_mt : SM_pt;}
 162
 163   /// indicates whether the split-merge orders on transverse mass or not.
 164   /// retained for backwards compatibility with 2.1.0b3
 165   bool split_merge_use_pt_weighted_splitting() const {return _use_pt_weighted_splitting;}
 166   void set_split_merge_use_pt_weighted_splitting(bool val) {
 167     _use_pt_weighted_splitting = val;}
 168
 169   // the things that are required by base class
 170   virtual std::string description () const;
 171   virtual void run_clustering(ClusterSequence &) const ;
 172
 173 protected:
 174   virtual void reset_stored_plugin() const;
 175
 176 private:
 177   double _protojet_ptmin;
 178   SplitMergeScale _split_merge_scale;
 179
 180   bool _use_pt_weighted_splitting;
 181
 182   // part needed for the cache
 183   // variables for caching the results and the input
 184   static std::auto_ptr<SISConePlugin          > stored_plugin;
 185   static std::auto_ptr<std::vector<PseudoJet> > stored_particles;
 186   static std::auto_ptr<siscone::Csiscone      > stored_siscone;
 187 };
 188
 189
 190 //======================================================================
 191 /// Class that provides extra information about a SISCone clustering
 192 class SISConeExtras : public SISConeBaseExtras {
 193 public:
 194   /// constructor
 195   //  it just initialises the pass information
 196   SISConeExtras(int nparticles)
 197     : SISConeBaseExtras(nparticles){}
 198
 199   /// access to the siscone jet def plugin (more convenient than
 200   /// getting it from the original jet definition, because here it's
 201   /// directly of the right type (rather than the base type)
 202   const SISConePlugin* jet_def_plugin() const {
 203     return dynamic_cast<const SISConePlugin*>(_jet_def_plugin);
 204   }
 205
 206 private:
 207   // let us be written to by SISConePlugin
 208   friend class SISConePlugin;
 209 };
 210
 211 } // fastjet namespace
 212
 213 #endif // __SISCONEPLUGIN_HH__
 214