![]() |
|
|||
File indexing completed on 2025-07-18 09:06:09
0001 //FJSTARTHEADER 0002 // $Id: PxConePlugin.hh 3433 2014-07-23 08:17:03Z salam $ 0003 // 0004 // Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez 0005 // 0006 //---------------------------------------------------------------------- 0007 // This file is part of FastJet. 0008 // 0009 // FastJet is free software; you can redistribute it and/or modify 0010 // it under the terms of the GNU General Public License as published by 0011 // the Free Software Foundation; either version 2 of the License, or 0012 // (at your option) any later version. 0013 // 0014 // The algorithms that underlie FastJet have required considerable 0015 // development. They are described in the original FastJet paper, 0016 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use 0017 // FastJet as part of work towards a scientific publication, please 0018 // quote the version you use and include a citation to the manual and 0019 // optionally also to hep-ph/0512210. 0020 // 0021 // FastJet is distributed in the hope that it will be useful, 0022 // but WITHOUT ANY WARRANTY; without even the implied warranty of 0023 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0024 // GNU General Public License for more details. 0025 // 0026 // You should have received a copy of the GNU General Public License 0027 // along with FastJet. If not, see <http://www.gnu.org/licenses/>. 0028 //---------------------------------------------------------------------- 0029 //FJENDHEADER 0030 0031 #ifndef __PXCONEPLUGIN_HH__ 0032 #define __PXCONEPLUGIN_HH__ 0033 0034 #include "fastjet/JetDefinition.hh" 0035 0036 // questionable whether this should be in fastjet namespace or not... 0037 0038 //FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 0039 namespace Rivet { 0040 //---------------------------------------------------------------------- 0041 // 0042 /// @ingroup plugins 0043 /// \class PxConePlugin 0044 /// Implementation of the PxCone algorithm (plugin for fastjet v2.1 upwards) 0045 /// 0046 /// PxConePlugin is a plugin for fastjet (v2.1 upwards) that provides 0047 /// an interface to the fortran pxcone iterative cone algorithm with 0048 /// midpoint seeds. 0049 /// 0050 /// Pxcone was written by Luis del Pozo and Michael H. Seymour. It is 0051 /// not a "supported" program, so if you encounter problems, you are 0052 /// on your own... 0053 /// 0054 /// Note that pxcone sometimes encounters non-stable iterations; in 0055 /// such cases it returns an error -- the plugin propagates this by 0056 /// throwing a fastjet::Error exception; if the user wishes to have 0057 /// robust code, they should catch this exception. 0058 /// 0059 /// Pxcone has a hard-coded limit (by default 4000) on the maximum 0060 /// number of particles and protojets; if the number of particles or 0061 /// protojets exceeds this, again a fastjet::Error exception will be 0062 /// thrown. 0063 /// 0064 /// The functionality of pxcone is described at 0065 /// http://www.hep.man.ac.uk/u/wplano/ConeJet.ps 0066 // 0067 //---------------------------------------------------------------------- 0068 class PxConePlugin : public fastjet::JetDefinition::Plugin { 0069 public: 0070 0071 /// constructor for the PxConePlugin, whose arguments have the 0072 /// following meaning: 0073 /// 0074 /// - the cone_radius is as usual in cone algorithms 0075 /// 0076 /// - stables cones (protojets) below min_jet_energy are discarded 0077 /// before calling the splitting procedure to resolve overlaps 0078 /// (called epslon in pxcone). 0079 /// 0080 /// - when two protojets overlap, if 0081 /// (overlapping_Et)/(Et_of_softer_protojet) < overlap_threshold 0082 /// the overlapping energy is split between the two protojets; 0083 /// otherwise the less energetic protojet is discarded. Called 0084 /// ovlim in pxcone. 0085 /// 0086 /// - pxcone carries out p-scheme recombination, and the resulting 0087 /// jets are massless; setting E_scheme_jets = true (default 0088 /// false) doesn't change the jet composition, but the final 0089 /// momentum sum for the jets is carried out by direct 0090 /// four-vector addition instead of p-scheme recombination. 0091 /// 0092 PxConePlugin (double cone_radius_in , 0093 double min_jet_energy_in = 5.0 , 0094 double overlap_threshold_in = 0.5, 0095 bool E_scheme_jets_in = false) : 0096 _cone_radius (cone_radius_in ), 0097 _min_jet_energy (min_jet_energy_in ), 0098 _overlap_threshold (overlap_threshold_in), 0099 _E_scheme_jets (E_scheme_jets_in ) { 0100 std::string msg = "Using own c++ version of PxCone, since FastJet doesn't install it by default. "; 0101 msg += "Please notify the Rivet authors if this behaviour should be changed."; 0102 std::cerr << msg << std::endl; 0103 } 0104 0105 0106 // some functions to return info about parameters ---------------- 0107 0108 /// the cone radius 0109 double cone_radius () const {return _cone_radius ;} 0110 0111 /// minimum jet energy (protojets below this are thrown own before 0112 /// merging/splitting) -- called epslon in pxcone 0113 double min_jet_energy () const {return _min_jet_energy ;} 0114 0115 /// Maximum fraction of overlap energy in a jet -- called ovlim in pxcone. 0116 double overlap_threshold () const {return _overlap_threshold ;} 0117 0118 /// if true then the final jets are returned as the E-scheme recombination 0119 /// of the particle momenta (by default, pxcone returns massless jets with 0120 /// a mean phi,eta type of recombination); regardless of what is 0121 /// returned, the internal pxcone jet-finding procedure is 0122 /// unaffected. 0123 bool E_scheme_jets() const {return _E_scheme_jets ;} 0124 0125 0126 // the things that are required by base class 0127 virtual std::string description () const; 0128 virtual void run_clustering(fastjet::ClusterSequence &) const; 0129 /// the plugin mechanism's standard way of accessing the jet radius 0130 virtual double R() const {return cone_radius();} 0131 0132 private: 0133 0134 double _cone_radius ; 0135 double _min_jet_energy ; 0136 double _overlap_threshold ; 0137 0138 bool _E_scheme_jets; 0139 0140 static bool _first_time; 0141 0142 /// print a banner for reference to the 3rd-party code 0143 void _print_banner(std::ostream *ostr) const; 0144 }; 0145 0146 0147 0148 // actual physical parameters: 0149 // 0150 // coner 0151 // epsilon 0152 // ovlim 0153 0154 void pxcone_ 0155 ( 0156 int mode , // 1=>e+e-, 2=>hadron-hadron 0157 int ntrak , // Number of particles 0158 int itkdm , // First dimension of PTRAK array: 0159 const double * ptrak , // Array of particle 4-momenta (Px,Py,Pz,E) 0160 double coner , // Cone size (half angle) in radians 0161 double epslon , // Minimum Jet energy (GeV) 0162 double ovlim , // Maximum fraction of overlap energy in a jet 0163 int mxjet , // Maximum possible number of jets 0164 int & njet , // Number of jets found 0165 double * pjet, // 5-vectors of jets 0166 int * ipass, // Particle k belongs to jet number IPASS(k)-1 0167 // IPASS = -1 if not assosciated to a jet 0168 int * ijmul, // Jet i contains IJMUL[i] particles 0169 int * ierr // = 0 if all is OK ; = -1 otherwise 0170 ); 0171 0172 //FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh 0173 } 0174 0175 #endif 0176
[ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |
![]() ![]() |