FastJet 3.0.2
|
00001 //STARTHEADER 00002 // $Id: EECambridgePlugin.cc 2577 2011-09-13 15:11:38Z salam $ 00003 // 00004 // Copyright (c) 2007-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 // fastjet stuff 00030 #include "fastjet/ClusterSequence.hh" 00031 #include "fastjet/EECambridgePlugin.hh" 00032 #include "fastjet/NNH.hh" 00033 00034 // other stuff 00035 #include <sstream> 00036 #include <limits> 00037 00038 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 00039 00040 using namespace std; 00041 00042 //---------------------------------------------------------------------- 00043 /// class to help run an e+e- Cambridge algorithm 00044 class EECamBriefJet { 00045 public: 00046 void init(const PseudoJet & jet) { 00047 double norm = 1.0/sqrt(jet.modp2()); 00048 nx = jet.px() * norm; 00049 ny = jet.py() * norm; 00050 nz = jet.pz() * norm; 00051 } 00052 00053 double distance(const EECamBriefJet * jet) const { 00054 double dij = 1 - nx*jet->nx 00055 - ny*jet->ny 00056 - nz*jet->nz; 00057 return dij; 00058 } 00059 00060 double beam_distance() const { 00061 return numeric_limits<double>::max(); 00062 } 00063 00064 private: 00065 double nx, ny, nz; 00066 }; 00067 00068 00069 string EECambridgePlugin::description () const { 00070 ostringstream desc; 00071 desc << "EECambridge plugin with ycut = " << ycut() ; 00072 return desc.str(); 00073 } 00074 00075 void EECambridgePlugin::run_clustering(ClusterSequence & cs) const { 00076 int njets = cs.jets().size(); 00077 NNH<EECamBriefJet> nnh(cs.jets()); 00078 00079 double Q2 = cs.Q2(); 00080 00081 while (njets > 0) { 00082 int i, j, k; 00083 // here we get a minimum based on the purely angular variable from the NNH class 00084 // (called dij there, but vij in the Cambridge article (which uses dij for 00085 // a kt distance...) 00086 double vij = nnh.dij_min(i, j); // i,j are return values... 00087 00088 // next we work out the dij (ee kt distance), and based on its 00089 // value decide whether we have soft-freezing (represented here by 00090 // a "Beam" clustering) or not 00091 double dij; 00092 if (j >= 0) { 00093 double scale = min(cs.jets()[i].E(), cs.jets()[j].E()); 00094 dij = 2 * vij * scale * scale; 00095 if (dij > Q2 * ycut()) { 00096 // we'll call the softer partner a "beam" jet 00097 if (cs.jets()[i].E() > cs.jets()[j].E()) std::swap(i,j); 00098 j = -1; 00099 } 00100 } else { 00101 // for the last particle left, just use yij = 1 00102 dij = Q2; 00103 } 00104 00105 if (j >= 0) { 00106 cs.plugin_record_ij_recombination(i, j, dij, k); 00107 nnh.merge_jets(i, j, cs.jets()[k], k); 00108 } else { 00109 cs.plugin_record_iB_recombination(i, dij); 00110 nnh.remove_jet(i); 00111 } 00112 njets--; 00113 } 00114 00115 } 00116 00117 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh