[u/mrichter/AliRoot.git] / JETAN / fastjet / fastjet / SISConePlugin.hh

#ifndef __SISCONEPLUGIN_HH__
#define __SISCONEPLUGIN_HH__

#include "SISConeBasePlugin.hh"

// forward declaration of the siscone classes we'll need
namespace siscone{
  class Csiscone;
}


namespace fastjet {      // defined in fastjet/internal/base.hh

//----------------------------------------------------------------------
//
/// SISConePlugin is a plugin for fastjet (v2.1 upwards) that provides
/// an interface to the seedless infrared safe cone jet finder by
/// Gregory Soyez and Gavin Salam.
///
/// SISCone uses geometrical techniques to exhaustively consider all
/// possible distinct cones. It then finds out which ones are stable
/// and sends the result to the Tevatron Run-II type split-merge
/// procedure for overlapping cones.
///
/// Four parameters govern the "physics" of the algorithm:
///
///  - the cone_radius (this should be self-explanatory!)
///
///  - the overlap_threshold is the parameter which dictates how much
///    two jets must overlap (pt_overlap/min(pt1,pt2)) if they are to be 
///    merged
///
///  - Not all particles are in stable cones in the first round of
///    searching for stable cones; one can therefore optionally have the
///    the jet finder carry out additional passes of searching for
///    stable cones among particles that were in no stable cone in
///    previous passes --- the maximum number of passes carried out is
///    n_pass_max. If this is zero then additional passes are carried
///    out until no new stable cones are found.
///
///  - Protojet ptmin: protojets that are below this ptmin
///    (default = 0) are discarded before each iteration of the
///    split-merge loop.
///
/// One parameter governs some internal algorithmic shortcuts: 
///
/// - if "caching" is turned on then the last event clustered by
///   siscone is stored -- if the current event is identical and the
///   cone_radius and n_pass_mass are identical, then the only part of
///   the clustering that needs to be rerun is the split-merge part,
///   leading to significant speed gains; there is a small (O(N) storage
///   and speed) penalty for caching, so it should be kept off
///   (default) if only a single overlap_threshold is used.
///
/// The final jets can be accessed by requestion the
/// inclusive_jets(...) from the ClusterSequence object. Note that
/// these PseudoJets have their user_index() set to the index of the
/// pass in which they were found (first pass = 0). NB: This does not
/// currently work for jets that consist of a single particle.
///
/// For further information on the details of the algorithm see the
/// SISCone paper, arXiv:0704.0292 [JHEP 0705:086,2007].
///
/// For documentation about the implementation, see the
/// siscone/doc/html/index.html file.
//
class SISConePlugin : public SISConeBasePlugin{
public:

  /// enum for the different split-merge scale choices;
  /// Note that order _must_ be the same as in siscone
  enum SplitMergeScale {SM_pt,     ///< transverse momentum (E-scheme), IR unsafe
                        SM_Et,     ///< transverse energy (E-scheme), not long. boost invariant
                                   ///< original run-II choice [may not be implemented]
                        SM_mt,     ///< transverse mass (E-scheme), IR safe except
                                   ///< in decays of two identical narrow heavy particles
                        SM_pttilde ///< pt-scheme pt = \sum_{i in jet} |p_{ti}|, should
                                   ///< be IR safe in all cases
  };


  /// Main constructor for the SISCone Plugin class.  
  ///
  /// Note: wrt version prior to 2.4 this constructor differs in that a 
  /// the default value has been removed for overlap_threshold. The
  /// former has been removed because the old default of 0.5 was found
  /// to be unsuitable in high-noise environments; so the user should
  /// now explicitly think about the value for this -- we recommend
  /// 0.75.
  ///
  SISConePlugin (double cone_radius,
                 double overlap_threshold,
                 int    n_pass_max = 0,
                 double protojet_ptmin = 0.0, 
                 bool   caching = false,
                 SplitMergeScale  split_merge_scale = SM_pttilde,
                 double split_merge_stopping_scale = 0.0){
    _cone_radius           = cone_radius;
    _overlap_threshold     = overlap_threshold;
    _n_pass_max            = n_pass_max;
    _protojet_ptmin        = protojet_ptmin;
    _caching               = caching;   
    _split_merge_scale     = split_merge_scale;
    _split_merge_stopping_scale = split_merge_stopping_scale;
    _ghost_sep_scale       = 0.0;
    _use_pt_weighted_splitting = false;}


  /// Backwards compatible constructor for the SISCone Plugin class
  SISConePlugin (double cone_radius,
                 double overlap_threshold,
                 int    n_pass_max,
                 double protojet_ptmin, 
                 bool   caching ,
                 bool   split_merge_on_transverse_mass){
    _cone_radius           = cone_radius;
    _overlap_threshold     = overlap_threshold;
    _n_pass_max            = n_pass_max;
    _protojet_ptmin        = protojet_ptmin;
    _caching               = caching;   
    _split_merge_stopping_scale = 0.0;
    _split_merge_scale     = split_merge_on_transverse_mass ? SM_mt : SM_pttilde;
    _ghost_sep_scale       = 0.0;}
  
  /// backwards compatible constructor for the SISCone Plugin class
  /// (avoid using this in future).
  SISConePlugin (double cone_radius,
                 double overlap_threshold,
                 int    n_pass_max,
                 bool   caching ) {
    _cone_radius           = cone_radius;
    _overlap_threshold     = overlap_threshold;
    _n_pass_max            = n_pass_max;
    _protojet_ptmin        = 0.0;
    _caching               = caching;   
    _split_merge_scale     = SM_mt;
    _split_merge_stopping_scale = 0.0;
    _ghost_sep_scale       = 0.0;
    _use_pt_weighted_splitting = false;}

  /// minimum pt for a protojet to be considered in the split-merge step
  /// of the algorithm
  double protojet_ptmin  () const {return _protojet_ptmin  ;}

  /// return the scale to be passed to SISCone as the protojet_ptmin
  /// -- if we have a ghost separation scale that is above the
  /// protojet_ptmin, then the ghost_separation_scale becomes the
  /// relevant one to use here
  double protojet_or_ghost_ptmin  () const {return std::max(_protojet_ptmin,
                                                            _ghost_sep_scale);}

  /// indicates scale used in split-merge
  SplitMergeScale split_merge_scale() const {return _split_merge_scale;}
  /// sets scale used in split-merge
  void set_split_merge_scale(SplitMergeScale sms) {_split_merge_scale = sms;}

  /// indicates whether the split-merge orders on transverse mass or not.
  /// retained for backwards compatibility with 2.1.0b3
  bool split_merge_on_transverse_mass() const {return _split_merge_scale == SM_mt ;}
  void set_split_merge_on_transverse_mass(bool val) {
    _split_merge_scale = val  ? SM_mt : SM_pt;}

  /// indicates whether the split-merge orders on transverse mass or not.
  /// retained for backwards compatibility with 2.1.0b3
  bool split_merge_use_pt_weighted_splitting() const {return _use_pt_weighted_splitting;}
  void set_split_merge_use_pt_weighted_splitting(bool val) {
    _use_pt_weighted_splitting = val;}

  // the things that are required by base class
  virtual std::string description () const;
  virtual void run_clustering(ClusterSequence &) const ;

protected:
  virtual void reset_stored_plugin() const;

private:
  double _protojet_ptmin;
  SplitMergeScale _split_merge_scale;

  bool _use_pt_weighted_splitting;

  // part needed for the cache 
  // variables for caching the results and the input
  static std::auto_ptr<SISConePlugin          > stored_plugin;
  static std::auto_ptr<std::vector<PseudoJet> > stored_particles;
  static std::auto_ptr<siscone::Csiscone      > stored_siscone;
};


//======================================================================
/// Class that provides extra information about a SISCone clustering
class SISConeExtras : public SISConeBaseExtras {
public:
  /// constructor
  //  it just initialises the pass information 
  SISConeExtras(int nparticles)
    : SISConeBaseExtras(nparticles){}

  /// access to the siscone jet def plugin (more convenient than
  /// getting it from the original jet definition, because here it's
  /// directly of the right type (rather than the base type)
  const SISConePlugin* jet_def_plugin() const {
    return dynamic_cast<const SISConePlugin*>(_jet_def_plugin);
  }

private:
  // let us be written to by SISConePlugin
  friend class SISConePlugin;
};

} // fastjet namespace 

#endif // __SISCONEPLUGIN_HH__
Commit	Line	Data
370be031	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