FastJet 3.0.2
|
00001 #ifndef __FASTJET_TOOLS_PRUNER_HH__ 00002 #define __FASTJET_TOOLS_PRUNER_HH__ 00003 00004 //STARTHEADER 00005 // $Id: Pruner.hh 2616 2011-09-30 18:03:40Z salam $ 00006 // 00007 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez 00008 // 00009 //---------------------------------------------------------------------- 00010 // This file is part of FastJet. 00011 // 00012 // FastJet is free software; you can redistribute it and/or modify 00013 // it under the terms of the GNU General Public License as published by 00014 // the Free Software Foundation; either version 2 of the License, or 00015 // (at your option) any later version. 00016 // 00017 // The algorithms that underlie FastJet have required considerable 00018 // development and are described in hep-ph/0512210. If you use 00019 // FastJet as part of work towards a scientific publication, please 00020 // include a citation to the FastJet paper. 00021 // 00022 // FastJet is distributed in the hope that it will be useful, 00023 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00024 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00025 // GNU General Public License for more details. 00026 // 00027 // You should have received a copy of the GNU General Public License 00028 // along with FastJet. If not, see <http://www.gnu.org/licenses/>. 00029 //---------------------------------------------------------------------- 00030 //ENDHEADER 00031 00032 #include "fastjet/ClusterSequence.hh" 00033 #include "fastjet/WrappedStructure.hh" 00034 #include "fastjet/tools/Transformer.hh" 00035 #include <iostream> 00036 #include <string> 00037 00038 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 00039 00040 // fwd declarations 00041 class Pruner; 00042 class PrunerStructure; 00043 class PruningRecombiner; 00044 class PruningPlugin; 00045 00046 //---------------------------------------------------------------------- 00047 /// @ingroup tools_generic 00048 /// \class Pruner 00049 /// Transformer that prunes a jet 00050 /// 00051 /// This transformer prunes a jet according to the ideas presented in 00052 /// arXiv:0903.5081 (S.D. Ellis, C.K. Vermilion and J.R. Walsh). 00053 /// 00054 /// The jet's constituents are reclustered with a user-specified jet 00055 /// definition, with the modification that objects i and j are only 00056 /// recombined if at least one of the following two criteria is 00057 /// satisfied: 00058 /// 00059 /// - the geometric distance between i and j is smaller than 'Rcut' 00060 /// with Rcut = Rcut_factor*2m/pt (Rcut_factor is a parameter of 00061 /// the Pruner and m and pt obtained from the jet being pruned) 00062 /// - the transverse momenta of i and j are at least 'zcut' p_t(i+j) 00063 /// 00064 /// If both these criteria fail, i and j are not recombined, the 00065 /// harder of i and j is kept, and the softer is rejected. 00066 /// 00067 /// Usage: 00068 /// \code 00069 /// Pruner pruner(jet_def, zcut, Rcut_factor); 00070 /// PseudoJet pruned_jet = pruner(jet); 00071 /// \endcode 00072 /// 00073 /// The pruned_jet has a valid associated cluster sequence. In addition 00074 /// the subjets of the original jet that have been vetoed by pruning 00075 /// (i.e. have been 'pruned away') can be accessed using 00076 /// 00077 /// \code 00078 /// vector<PseudoJet> rejected_subjets = pruned_jet.structure_of<Pruner>().rejected(); 00079 /// \endcode 00080 /// 00081 /// If the re-clustering happens to find more than a single inclusive 00082 /// jet (this should normally not happen if the radius of the jet 00083 /// definition used for the reclustering was set large enough), 00084 /// the hardest of these jets is retured as the result of the 00085 /// Pruner. The other jets can be accessed through 00086 /// 00087 /// \code 00088 /// vector<PseudoJet> extra_jets = pruned_jet.structure_of<Pruner>().extra_jets(); 00089 /// \endcode 00090 /// 00091 /// Instead of using Rcut_factor and zcut, one can alternatively 00092 /// construct a Pruner by passing two (pointers to) functions of 00093 /// PseudoJet that dynamically compute the Rcut and zcut to 00094 /// be used for the jet being pruned. 00095 /// 00096 /// When the jet being pruned has area support and explicit ghosts, 00097 /// the resulting pruned jet will likewise have area. 00098 /// 00099 //---------------------------------------------------------------------- 00100 class Pruner : public Transformer{ 00101 public: 00102 /// minimal constructor, which takes a jet algorithm, sets the radius 00103 /// to JetDefinition::max_allowable_R (practically equivalent to 00104 /// infinity) and also tries to use a recombiner based on the one in 00105 /// the jet definition of the particular jet being pruned. 00106 /// 00107 /// \param jet_alg the jet algorithm for the internal clustering 00108 /// \param zcut pt-fraction cut in the pruning 00109 /// \param Rcut_factor the angular distance cut in the pruning will be 00110 /// Rcut_factor * 2m/pt 00111 Pruner(const JetAlgorithm jet_alg, double zcut, double Rcut_factor) 00112 : _jet_def(jet_alg, JetDefinition::max_allowable_R), 00113 _zcut(zcut), _Rcut_factor(Rcut_factor), 00114 _zcut_dyn(0), _Rcut_dyn(0), _get_recombiner_from_jet(true) {} 00115 00116 00117 /// alternative ctor in which the full reclustering jet definition can 00118 /// be specified. 00119 /// 00120 /// \param jet_def the jet definition for the internal clustering 00121 /// \param zcut pt-fraction cut in the pruning 00122 /// \param Rcut_factor the angular distance cut in the pruning will be 00123 /// Rcut_factor * 2m/pt 00124 Pruner(const JetDefinition &jet_def, double zcut, double Rcut_factor) 00125 : _jet_def(jet_def), 00126 _zcut(zcut), _Rcut_factor(Rcut_factor), 00127 _zcut_dyn(0), _Rcut_dyn(0), _get_recombiner_from_jet(false) {} 00128 00129 00130 /// alternative ctor in which the pt-fraction cut and angular distance 00131 /// cut are functions of the jet being pruned. 00132 /// 00133 /// \param jet_def the jet definition for the internal clustering 00134 /// \param zcut_dyn dynamic pt-fraction cut in the pruning 00135 /// \param Rcut_dyn dynamic angular distance cut in the pruning 00136 Pruner(const JetDefinition &jet_def, 00137 FunctionOfPseudoJet<double> *zcut_dyn, 00138 FunctionOfPseudoJet<double> *Rcut_dyn); 00139 00140 /// action on a single jet 00141 virtual PseudoJet result(const PseudoJet &jet) const; 00142 00143 /// description 00144 virtual std::string description() const; 00145 00146 // the type of the associated structure 00147 typedef PrunerStructure StructureType; 00148 00149 private: 00150 /// check if the jet has explicit_ghosts (knowing that there is an 00151 /// area support) 00152 bool _check_explicit_ghosts(const PseudoJet &jet) const; 00153 00154 /// return a pointer to a "common" recombiner if there is one, 00155 /// alternatively a null pointer. 00156 const JetDefinition::Recombiner * _get_common_recombiner(const PseudoJet &jet) const; 00157 00158 JetDefinition _jet_def; ///< the internal jet definition 00159 double _zcut; ///< the pt-fraction cut 00160 double _Rcut_factor; ///< the angular separation cut factor 00161 FunctionOfPseudoJet<double> *_zcut_dyn; ///< dynamic zcut 00162 FunctionOfPseudoJet<double> *_Rcut_dyn; ///< dynamic Rcut 00163 bool _get_recombiner_from_jet; ///< true for minimal constructor, 00164 ///< causes recombiner to be set equal 00165 ///< to that already used in the jet 00166 ///< (if it can be deduced) 00167 }; 00168 00169 00170 //---------------------------------------------------------------------- 00171 /// @ingroup tools_generic 00172 /// \class PrunerStructure 00173 /// The structure associated with a PseudoJet thas has gone through a 00174 /// Pruner transformer 00175 //---------------------------------------------------------------------- 00176 class PrunerStructure : public WrappedStructure{ 00177 public: 00178 /// default ctor 00179 /// \param result_jet the jet for which we have to keep the structure 00180 PrunerStructure(const PseudoJet & result_jet) 00181 : WrappedStructure(result_jet.structure_shared_ptr()){} 00182 00183 /// description 00184 virtual std::string description() const{ return "Pruned PseudoJet";} 00185 00186 /// return the constituents that have been rejected 00187 std::vector<PseudoJet> rejected() const{ 00188 return validated_cs()->childless_pseudojets(); 00189 } 00190 00191 /// return the other jets that may have been found along with the 00192 /// result of the pruning 00193 /// The resulting vector is sorted in pt 00194 std::vector<PseudoJet> extra_jets() const; 00195 00196 protected: 00197 friend class Pruner; ///< to allow setting the internal information 00198 }; 00199 00200 //---------------------------------------------------------------------- 00201 /// \if internal_doc 00202 /// @ingroup internal 00203 /// \class PruningRecombiner 00204 /// recombines the objects that are not vetoed by pruning 00205 /// 00206 /// This recombiner only recombines, using the provided 'recombiner', 00207 /// objects (i and j) that pass at least one of the following two criteria: 00208 /// 00209 /// - the geometric distance between i and j is smaller than 'Rcut' 00210 /// - the transverse momenta of i and j are at least 'zcut' p_t(i+j) 00211 /// 00212 /// If both these criteria fail, the hardest of i and j is kept and 00213 /// the softest is rejected. 00214 /// 00215 /// Note that this in not meant for standalone use [in particular 00216 /// because it could lead to memory issues due to the rejected indices 00217 /// stored internally]. 00218 /// 00219 /// \endif 00220 class PruningRecombiner : public JetDefinition::Recombiner{ 00221 public: 00222 /// ctor 00223 /// \param zcut transverse momentum fraction cut 00224 /// \param Rcut angular separation cut 00225 /// \param recomb pointer to a recombiner to use to cluster pairs 00226 PruningRecombiner(double zcut, double Rcut, 00227 const JetDefinition::Recombiner *recombiner) 00228 : _zcut2(zcut*zcut), _Rcut2(Rcut*Rcut), 00229 _recombiner(recombiner){} 00230 00231 /// perform a recombination taking into account the pruning 00232 /// conditions 00233 virtual void recombine(const PseudoJet &pa, 00234 const PseudoJet &pb, 00235 PseudoJet &pab) const; 00236 00237 /// returns the description of the recombiner 00238 virtual std::string description() const; 00239 00240 /// return the history indices that have been pruned away 00241 const std::vector<unsigned int> & rejected() const{ return _rejected;} 00242 00243 /// clears the list of rejected indices 00244 /// 00245 /// If one decides to use this recombiner standalone, one has to 00246 /// call this after each clustering in order for the rejected() vector 00247 /// to remain sensible and not grow to infinite size. 00248 void clear_rejected(){ _rejected.clear();} 00249 00250 private: 00251 double _zcut2; ///< transverse momentum fraction cut (squared) 00252 double _Rcut2; ///< angular separation cut (squared) 00253 const JetDefinition::Recombiner *_recombiner; ///< the underlying recombiner to use 00254 mutable std::vector<unsigned int> _rejected; ///< list of rejected history indices 00255 }; 00256 00257 00258 //---------------------------------------------------------------------- 00259 /// \if internal_doc 00260 /// @ingroup internal 00261 /// \class PruningPlugin 00262 /// FastJet internal plugin that clusters the particles using the 00263 /// PruningRecombiner. 00264 /// 00265 /// See PruningRecombiner for a description of what pruning does. 00266 /// 00267 /// Note that this is an internal FastJet class used by the Pruner 00268 /// transformer and it is not meant to be used as a standalone clustering 00269 /// tool. 00270 /// 00271 /// \endif 00272 //---------------------------------------------------------------------- 00273 class PruningPlugin : public JetDefinition::Plugin{ 00274 public: 00275 /// ctor 00276 /// \param jet_def the jet definition to be used for the 00277 /// internal clustering 00278 /// \param zcut transverse momentum fraction cut 00279 /// \param Rcut angular separation cut 00280 PruningPlugin(const JetDefinition &jet_def, double zcut, double Rcut) 00281 : _jet_def(jet_def), _zcut(zcut), _Rcut(Rcut){} 00282 00283 /// the actual clustering work for the plugin 00284 virtual void run_clustering(ClusterSequence &input_cs) const; 00285 00286 /// description of the plugin 00287 virtual std::string description() const; 00288 00289 /// returns the radius 00290 virtual double R() const {return _jet_def.R();} 00291 00292 private: 00293 JetDefinition _jet_def; ///< the internal jet definition 00294 double _zcut; ///< transverse momentum fraction cut 00295 double _Rcut; ///< angular separation cut 00296 }; 00297 00298 00299 00300 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh 00301 00302 #endif // __FASTJET_TOOLS_PRUNER_HH__