FastJet  3.3.1
ATLASConePlugin.cc
1 //FJSTARTHEADER
2 // $Id: ATLASConePlugin.cc 4354 2018-04-22 07:12:37Z salam $
3 //
4 // Copyright (c) 2007-2018, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5 //
6 //----------------------------------------------------------------------
7 // This file is part of FastJet.
8 //
9 // FastJet is free software; you can redistribute it and/or modify
10 // it under the terms of the GNU General Public License as published by
11 // the Free Software Foundation; either version 2 of the License, or
12 // (at your option) any later version.
13 //
14 // The algorithms that underlie FastJet have required considerable
15 // development. They are described in the original FastJet paper,
16 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17 // FastJet as part of work towards a scientific publication, please
18 // quote the version you use and include a citation to the manual and
19 // optionally also to hep-ph/0512210.
20 //
21 // FastJet is distributed in the hope that it will be useful,
22 // but WITHOUT ANY WARRANTY; without even the implied warranty of
23 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 // GNU General Public License for more details.
25 //
26 // You should have received a copy of the GNU General Public License
27 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28 //----------------------------------------------------------------------
29 //FJENDHEADER
30 
31 // fastjet stuff
32 #include "fastjet/ClusterSequence.hh"
33 #include "fastjet/ATLASConePlugin.hh"
34 
35 // SpartyJet stuff
36 #include "CommonUtils.hh"
37 #include "JetConeFinderTool.hh"
38 #include "JetSplitMergeTool.hh"
39 
40 // other stuff
41 #include <vector>
42 #include <sstream>
43 
44 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
45 
46 using namespace std;
47 
48 bool ATLASConePlugin::_first_time = true;
49 
50 string ATLASConePlugin::description () const {
51  ostringstream desc;
52  desc << "ATLASCone plugin with R = "<< _radius
53  << ", seed threshold = " << _seedPt
54  << ", overlap threshold f = " << _f;
55  return desc.str();
56 }
57 
58 void ATLASConePlugin::run_clustering(ClusterSequence & clust_seq) const {
59  // print a banner if we run this for the first time
60  _print_banner(clust_seq.fastjet_banner_stream());
61 
62 
63  // transfer the list of PseudoJet into a atlas::Jet::jet_list_t
64  // jet_list_t is a vector<Jet*>
65  // We set the index of the 4-vect to trace the constituents at the end
66  //------------------------------------------------------------------
67  // cout << "ATLASConePlugin: transferring vectors from ClusterSequence" << endl;
68  atlas::JetConeFinderTool::jetcollection_t jets_ptr;
69  vector<atlas::Jet*> particles_ptr;
70 
71  for (unsigned int i=0 ; i<clust_seq.jets().size() ; i++) {
72  const PseudoJet & mom = clust_seq.jets()[i];
73 
74  // first create the particle
75  atlas::Jet *particle = new atlas::Jet(mom.px(), mom.py(), mom.pz(), mom.E(), i);
76  particles_ptr.push_back(particle);
77 
78  // then add it to the list of particles we'll use for the clustering
79  atlas::Jet *jet = new atlas::Jet;
80  jet->set_index(particle->index());
81  jet->addConstituent(particle);
82 
83  // and finally add that one to the list of jets
84  jets_ptr.push_back(jet);
85  }
86  // cout << "ATLASCone: " << jets_ptr.size() << " particles to cluster" << endl;
87 
88  // search the stable cones
89  //------------------------------------------------------------------
90  // cout << "ATLASConePlugin: searching for stable cones" << endl;
91  atlas::JetConeFinderTool stable_cone_finder;
92 
93  // set the parameters
94  stable_cone_finder.m_coneR = _radius;
95  stable_cone_finder.m_seedPt = _seedPt;
96 
97  // really do the search.
98  // Note that the list of protocones is returned
99  // through the argument
100  stable_cone_finder.execute(jets_ptr);
101  // cout << "ATLASCone: " << jets_ptr.size() << " stable cones found" << endl;
102 
103  // perform the split-merge
104  //------------------------------------------------------------------
105  // cout << "ATLASConePlugin: running the split-merge" << endl;
106  atlas::JetSplitMergeTool split_merge;
107 
108  // set the parameters
109  split_merge.m_f = _f;
110 
111  // do the work
112  // again, the list of jets is returned through the argument
113  split_merge.execute(&jets_ptr);
114  // cout << "ATLASCone: " << jets_ptr.size() << " jets after split--merge" << endl;
115 
116  // build the FastJet jets (a la SISConePlugin)
117  //------------------------------------------------------------------
118  // cout << "ATLASConePlugin: backporting jets to the ClusterSequence" << endl;
119  for (atlas::Jet::jet_list_t::iterator jet_it = jets_ptr.begin() ;
120  jet_it != jets_ptr.end(); jet_it++){
121  // iterate over the constituents, starting from the first one
122  // that we just take as a reference
123  atlas::Jet::constit_vect_t::iterator constit_it = (*jet_it)->firstConstituent();
124  // cout << " atlas: jet has " << (*jet_it)->getConstituentNum() << " constituents" << endl;
125  int jet_k = (*constit_it)->index();
126  constit_it++;
127 
128  // loop over the remaining particles
129  while (constit_it != (*jet_it)->lastConstituent()){
130  // take the last result of the merge
131  int jet_i = jet_k;
132  // and the next element of the jet
133  int jet_j = (*constit_it)->index();
134  // and merge them (with a fake dij)
135  double dij = 0.0;
136 
137  // create the new jet by hand so that we can adjust its user index
138  // Note again the use of the E-scheme recombination here!
139  PseudoJet newjet = clust_seq.jets()[jet_i] + clust_seq.jets()[jet_j];
140  clust_seq.plugin_record_ij_recombination(jet_i, jet_j, dij, newjet, jet_k);
141 
142  // jump to the next constituent
143  constit_it++;
144  }
145 
146  // we have merged all the jet's particles into a single object, so now
147  // "declare" it to be a beam (inclusive) jet.
148  // [NB: put a sensible looking d_iB just to be nice...]
149  double d_iB = clust_seq.jets()[jet_k].perp2();
150  clust_seq.plugin_record_iB_recombination(jet_k, d_iB);
151  }
152 
153  // cout << "ATLASConePlugin: Bye" << endl;
154  clear_list(particles_ptr);
155  clear_list(jets_ptr);
156 }
157 
158 // print a banner for reference to the 3rd-party code
159 void ATLASConePlugin::_print_banner(ostream *ostr) const{
160  if (! _first_time) return;
161  _first_time=false;
162 
163  // make sure the user has not set the banner stream to NULL
164  if (!ostr) return;
165 
166  (*ostr) << "#-------------------------------------------------------------------------" << endl;
167  (*ostr) << "# You are running the ATLAS Cone plugin for FastJet " << endl;
168  (*ostr) << "# Original code from SpartyJet; interface by the FastJet authors " << endl;
169  (*ostr) << "# If you use this plugin, please cite " << endl;
170  (*ostr) << "# P.A. Delsart, K. Geerlings, J. Huston, B. Martin and C. Vermilion, " << endl;
171  (*ostr) << "# SpartyJet, http://projects.hepforge.org/spartyjet " << endl;
172  (*ostr) << "# in addition to the usual FastJet reference. " << endl;
173  (*ostr) << "#-------------------------------------------------------------------------" << endl;
174 
175  // make sure we really have the output done.
176  ostr->flush();
177 }
178 
179 
180 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
deals with clustering
void plugin_record_ij_recombination(int jet_i, int jet_j, double dij, int &newjet_k)
record the fact that there has been a recombination between jets()[jet_i] and jets()[jet_k], with the specified dij, and return the index (newjet_k) allocated to the new jet, whose momentum is assumed to be the 4-vector sum of that of jet_i and jet_j
void plugin_record_iB_recombination(int jet_i, double diB)
record the fact that there has been a recombination between jets()[jet_i] and the beam...
static std::ostream * fastjet_banner_stream()
returns a pointer to the stream to be used to print banners (cout by default).
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
const std::vector< PseudoJet > & jets() const
allow the user to access the internally stored _jets() array, which contains both the initial particl...