FastJet  3.3.3
NestedDefsPlugin.cc
1 //FJSTARTHEADER
2 // $Id: NestedDefsPlugin.cc 4420 2019-11-29 09:28:20Z soyez $
3 //
4 // Copyright (c) 2007-2019, 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 // TODO
32 // ? Maybe one could provide additional recomb. dists as an "extra".;
33 
34 // fastjet stuff
35 #include "fastjet/ClusterSequence.hh"
36 #include "fastjet/NestedDefsPlugin.hh"
37 
38 // other stuff
39 #include <vector>
40 #include <sstream>
41 
42 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
43 
44 using namespace std;
45 
46 string NestedDefsPlugin::description () const {
47  ostringstream desc;
48 
49  desc << "NestedDefs: successive application of " ;
50  unsigned int i=1;
51  for (list<JetDefinition>::const_iterator it=_defs.begin();it!=_defs.end();it++){
52  desc << "Definition " << i++ << " [" << it->description() << "] - ";
53  }
54 
55  return desc.str();
56 }
57 
58 void NestedDefsPlugin::run_clustering(ClusterSequence & clust_seq) const {
59  vector<PseudoJet> momenta;
60 
61  // build the initial list of particles
62  momenta = clust_seq.jets();
63  unsigned int step_n = momenta.size();
64 
65  // initialise the conversion table, which works as follows
66  // conversion_table[step_cs_jet_index] = main_cs_jet_index
67  vector<unsigned int> conversion_table(2*step_n);
68  vector<unsigned int> new_conversion_table;
69  for (unsigned int i=0;i<step_n;i++)
70  conversion_table[i]=i;
71 
72  // Now the steps go as follows:
73  // for each definition in the list,
74  // - do the clustering,
75  // - copy the history into the main one
76  // - update the list of momenta and the index conversion table
77  list<JetDefinition>::const_iterator def_iterator = _defs.begin();
78  unsigned int def_index=0;
79  bool last_def=false;
80 
81  while (def_iterator!=_defs.end()){
82  last_def = (def_index == (_defs.size()-1));
83 
84  // do the clustering
85  ClusterSequence step_cs(momenta, *def_iterator);
86 
87  // clear the momenta as we shall fill them again
88  momenta.clear();
89  new_conversion_table.clear();
90 
91  // retrieve the history
92  const vector<ClusterSequence::history_element> & step_history = step_cs.history();
93 
94  // copy the history
95  // note that we skip the initial steps which are just the
96  // declaration of the particles.
97  vector<ClusterSequence::history_element>::const_iterator
98  hist_iterator = step_history.begin();
99 
100  for (unsigned int i=step_n;i!=0;i--)
101  hist_iterator++;
102 
103  while (hist_iterator != step_history.end()){
104  // check if it is a recombination with the beam or a simple recombination
105  if (hist_iterator->parent2 == ClusterSequence::BeamJet){
106  // save this jet for future clustering
107  // unless we've reached the last def in which case, record the clustering
108  unsigned int step_jet_index = step_cs.history()[hist_iterator->parent1].jetp_index;
109  if (last_def){
110  clust_seq.plugin_record_iB_recombination(conversion_table[step_jet_index],
111  hist_iterator->dij);
112  } else {
113  momenta.push_back(step_cs.jets()[step_jet_index]);
114  new_conversion_table.push_back(conversion_table[step_jet_index]);
115  }
116  } else {
117  // record combination
118  // note that we set the recombination distance to 0 except for the last alg
119  unsigned int step_jet1_index = step_cs.history()[hist_iterator->parent1].jetp_index;
120  unsigned int step_jet2_index = step_cs.history()[hist_iterator->parent2].jetp_index;
121  PseudoJet newjet = step_cs.jets()[hist_iterator->jetp_index];
122  int jet_k;
123  clust_seq.plugin_record_ij_recombination(conversion_table[step_jet1_index],
124  conversion_table[step_jet2_index],
125  last_def ? hist_iterator->dij : 0.0,
126  newjet, jet_k);
127 
128  // save info in the conversion table for tracking purposes
129  conversion_table[hist_iterator->jetp_index]=jet_k;
130  }
131 
132  // go to the next history element
133  hist_iterator++;
134  }
135 
136  // finalise this step:
137  // - update nr of particles
138  // - update conversion table
139  step_n = momenta.size();
140  for (unsigned int i=0;i<step_n;i++)
141  conversion_table[i] = new_conversion_table[i];
142 
143  // go to the next alg
144  def_index++;
145  def_iterator++;
146  }
147 
148 }
149 
150 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
deals with clustering
const std::vector< history_element > & history() const
allow the user to access the raw internal history.
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...
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...