FastJet 3.0.2
|
00001 //STARTHEADER 00002 // $Id: ATLASConePlugin.cc 2777 2011-11-25 15:01:56Z soyez $ 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/ATLASConePlugin.hh" 00032 00033 // SpartyJet stuff 00034 #include "CommonUtils.hh" 00035 #include "JetConeFinderTool.hh" 00036 #include "JetSplitMergeTool.hh" 00037 00038 // other stuff 00039 #include <vector> 00040 #include <sstream> 00041 00042 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 00043 00044 using namespace std; 00045 00046 bool ATLASConePlugin::_first_time = true; 00047 00048 string ATLASConePlugin::description () const { 00049 ostringstream desc; 00050 desc << "ATLASCone plugin with R = "<< _radius 00051 << ", seed threshold = " << _seedPt 00052 << ", overlap threshold f = " << _f; 00053 return desc.str(); 00054 } 00055 00056 void ATLASConePlugin::run_clustering(ClusterSequence & clust_seq) const { 00057 // print a banner if we run this for the first time 00058 _print_banner(clust_seq.fastjet_banner_stream()); 00059 00060 00061 // transfer the list of PseudoJet into a atlas::Jet::jet_list_t 00062 // jet_list_t is a vector<Jet*> 00063 // We set the index of the 4-vect to trace the constituents at the end 00064 //------------------------------------------------------------------ 00065 // cout << "ATLASConePlugin: transferring vectors from ClusterSequence" << endl; 00066 atlas::JetConeFinderTool::jetcollection_t jets_ptr; 00067 vector<atlas::Jet*> particles_ptr; 00068 00069 for (unsigned int i=0 ; i<clust_seq.jets().size() ; i++) { 00070 const PseudoJet & mom = clust_seq.jets()[i]; 00071 00072 // first create the particle 00073 atlas::Jet *particle = new atlas::Jet(mom.px(), mom.py(), mom.pz(), mom.E(), i); 00074 particles_ptr.push_back(particle); 00075 00076 // then add it to the list of particles we'll use for the clustering 00077 atlas::Jet *jet = new atlas::Jet; 00078 jet->set_index(particle->index()); 00079 jet->addConstituent(particle); 00080 00081 // and finally add that one to the list of jets 00082 jets_ptr.push_back(jet); 00083 } 00084 // cout << "ATLASCone: " << jets_ptr.size() << " particles to cluster" << endl; 00085 00086 // search the stable cones 00087 //------------------------------------------------------------------ 00088 // cout << "ATLASConePlugin: searching for stable cones" << endl; 00089 atlas::JetConeFinderTool stable_cone_finder; 00090 00091 // set the parameters 00092 stable_cone_finder.m_coneR = _radius; 00093 stable_cone_finder.m_seedPt = _seedPt; 00094 00095 // really do the search. 00096 // Note that the list of protocones is returned 00097 // through the argument 00098 stable_cone_finder.execute(jets_ptr); 00099 // cout << "ATLASCone: " << jets_ptr.size() << " stable cones found" << endl; 00100 00101 // perform the split-merge 00102 //------------------------------------------------------------------ 00103 // cout << "ATLASConePlugin: running the split-merge" << endl; 00104 atlas::JetSplitMergeTool split_merge; 00105 00106 // set the parameters 00107 split_merge.m_f = _f; 00108 00109 // do the work 00110 // again, the list of jets is returned through the argument 00111 split_merge.execute(&jets_ptr); 00112 // cout << "ATLASCone: " << jets_ptr.size() << " jets after split--merge" << endl; 00113 00114 // build the FastJet jets (a la SISConePlugin) 00115 //------------------------------------------------------------------ 00116 // cout << "ATLASConePlugin: backporting jets to the ClusterSequence" << endl; 00117 for (atlas::Jet::jet_list_t::iterator jet_it = jets_ptr.begin() ; 00118 jet_it != jets_ptr.end(); jet_it++){ 00119 // iterate over the constituents, starting from the first one 00120 // that we just take as a reference 00121 atlas::Jet::constit_vect_t::iterator constit_it = (*jet_it)->firstConstituent(); 00122 // cout << " atlas: jet has " << (*jet_it)->getConstituentNum() << " constituents" << endl; 00123 int jet_k = (*constit_it)->index(); 00124 constit_it++; 00125 00126 // loop over the remaining particles 00127 while (constit_it != (*jet_it)->lastConstituent()){ 00128 // take the last result of the merge 00129 int jet_i = jet_k; 00130 // and the next element of the jet 00131 int jet_j = (*constit_it)->index(); 00132 // and merge them (with a fake dij) 00133 double dij = 0.0; 00134 00135 // create the new jet by hand so that we can adjust its user index 00136 // Note again the use of the E-scheme recombination here! 00137 PseudoJet newjet = clust_seq.jets()[jet_i] + clust_seq.jets()[jet_j]; 00138 clust_seq.plugin_record_ij_recombination(jet_i, jet_j, dij, newjet, jet_k); 00139 00140 // jump to the next constituent 00141 constit_it++; 00142 } 00143 00144 // we have merged all the jet's particles into a single object, so now 00145 // "declare" it to be a beam (inclusive) jet. 00146 // [NB: put a sensible looking d_iB just to be nice...] 00147 double d_iB = clust_seq.jets()[jet_k].perp2(); 00148 clust_seq.plugin_record_iB_recombination(jet_k, d_iB); 00149 } 00150 00151 // cout << "ATLASConePlugin: Bye" << endl; 00152 clear_list(particles_ptr); 00153 clear_list(jets_ptr); 00154 } 00155 00156 // print a banner for reference to the 3rd-party code 00157 void ATLASConePlugin::_print_banner(ostream *ostr) const{ 00158 if (! _first_time) return; 00159 _first_time=false; 00160 00161 // make sure the user has not set the banner stream to NULL 00162 if (!ostr) return; 00163 00164 (*ostr) << "#-------------------------------------------------------------------------" << endl; 00165 (*ostr) << "# You are running the ATLAS Cone plugin for FastJet " << endl; 00166 (*ostr) << "# Original code from SpartyJet; interface by the FastJet authors " << endl; 00167 (*ostr) << "# If you use this plugin, please cite " << endl; 00168 (*ostr) << "# P.A. Delsart, K. Geerlings, J. Huston, B. Martin and C. Vermilion, " << endl; 00169 (*ostr) << "# SpartyJet, http://projects.hepforge.org/spartyjet " << endl; 00170 (*ostr) << "# in addition to the usual FastJet reference. " << endl; 00171 (*ostr) << "#-------------------------------------------------------------------------" << endl; 00172 00173 // make sure we really have the output done. 00174 ostr->flush(); 00175 } 00176 00177 00178 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh