1 //FJSTARTHEADER
2 // $Id: CDFMidPointPlugin.cc 4442 2020-05-05 07:50:11Z soyez $
3 //
4 // Copyright (c) 2005-2020, 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 #include "fastjet/CDFMidPointPlugin.hh"
32 #include "fastjet/ClusterSequence.hh"
33 #include "fastjet/Error.hh"
34 #include <sstream>
35
36 // CDF stuff
37 #include "MidPointAlgorithm.hh"
38 #include "PhysicsTower.hh"
39 #include "Cluster.hh"
40
41 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
42
43 using namespace std;
44 using namespace cdf;
45
46 bool CDFMidPointPlugin::_first_time = true;
47
description() const48 string CDFMidPointPlugin::description () const {
49 ostringstream desc;
50
51 string sm_scale_string = "split-merge uses ";
52 switch(_sm_scale) {
53 case SM_pt:
54 sm_scale_string += "pt";
55 break;
56 case SM_Et:
57 sm_scale_string += "Et";
58 break;
59 case SM_mt:
60 sm_scale_string += "mt";
61 break;
62 case SM_pttilde:
63 sm_scale_string += "pttilde (scalar sum of pts)";
64 break;
65 default:
66 ostringstream err;
67 err << "Unrecognized split-merge scale choice = " << _sm_scale;
68 throw Error(err.str());
69 }
70
71
72 if (cone_area_fraction() == 1) {
73 desc << "CDF MidPoint jet algorithm, with " ;
74 } else {
75 desc << "CDF MidPoint+Searchcone jet algorithm, with ";
76 }
77 desc << "seed_threshold = " << seed_threshold () << ", "
78 << "cone_radius = " << cone_radius () << ", "
79 << "cone_area_fraction = " << cone_area_fraction () << ", "
80 << "max_pair_size = " << max_pair_size () << ", "
81 << "max_iterations = " << max_iterations () << ", "
82 << "overlap_threshold = " << overlap_threshold () << ", "
83 << sm_scale_string ;
84
85 return desc.str();
86 }
87
88
run_clustering(ClusterSequence & clust_seq) const89 void CDFMidPointPlugin::run_clustering(ClusterSequence & clust_seq) const {
90 // print a banner if we run this for the first time
91 _print_banner(clust_seq.fastjet_banner_stream());
92
93 // create the physics towers needed by the CDF code
94 vector<PhysicsTower> towers;
95 towers.reserve(clust_seq.jets().size());
96 for (unsigned i = 0; i < clust_seq.jets().size(); i++) {
97 LorentzVector fourvect(clust_seq.jets()[i].px(),
98 clust_seq.jets()[i].py(),
99 clust_seq.jets()[i].pz(),
100 clust_seq.jets()[i].E());
101 PhysicsTower tower(fourvect);
102 // misuse one of the indices for tracking, since the MidPoint
103 // implementation doesn't seem to make use of these indices
104 tower.calTower.iEta = i;
105 towers.push_back(tower);
106 }
107
108 // prepare the CDF algorithm
109 MidPointAlgorithm m(_seed_threshold,_cone_radius,_cone_area_fraction,
110 _max_pair_size,_max_iterations,_overlap_threshold,
111 MidPointAlgorithm::SplitMergeScale(_sm_scale));
112
113 // run the CDF algorithm
114 std::vector<Cluster> jets;
115 m.run(towers,jets);
116
117
118 // now transfer the jets back into our own structure -- we will
119 // mimic the cone code with a sequential recombination sequence in
120 // which the jets are built up by adding one particle at a time
121 for(vector<Cluster>::const_iterator jetIter = jets.begin();
122 jetIter != jets.end(); jetIter++) {
123 const vector<PhysicsTower> & tower_list = jetIter->towerList;
124 int jet_k = tower_list[0].calTower.iEta;
125
126 int ntow = int(jetIter->towerList.size());
127 for (int itow = 1; itow < ntow; itow++) {
128 int jet_i = jet_k;
129 // retrieve our misappropriated index for the jet
130 int jet_j = tower_list[itow].calTower.iEta;
131 // do a fake recombination step with dij=0
132 double dij = 0.0;
133 clust_seq.plugin_record_ij_recombination(jet_i, jet_j, dij, jet_k);
134 }
135
136 // NB: put a sensible looking d_iB just to be nice...
137 double d_iB = clust_seq.jets()[jet_k].perp2();
138 clust_seq.plugin_record_iB_recombination(jet_k, d_iB);
139 }
140
141
142 // following code is for testing only
143 //cout << endl;
144 //for(vector<Cluster>::const_iterator jetIter = jets.begin();
145 // jetIter != jets.end(); jetIter++) {
146 // cout << jetIter->fourVector.pt() << " " << jetIter->fourVector.y() << endl;
147 //}
148 //cout << "-----------------------------------------------------\n";
149 //vector<PseudoJet> ourjets(clust_seq.inclusive_jets());
150 //for (vector<PseudoJet>::const_reverse_iterator ourjet = ourjets.rbegin();
151 // ourjet != ourjets.rend(); ourjet++) {
152 // cout << ourjet->perp() << " " << ourjet->rap() << endl;
153 //}
154 //cout << endl;
155 }
156
157 // print a banner for reference to the 3rd-party code
_print_banner(ostream * ostr) const158 void CDFMidPointPlugin::_print_banner(ostream *ostr) const{
159 if (! _first_time) return;
160 _first_time=false;
161
162 // make sure the user has not set the banner stream to NULL
163 if (!ostr) return;
164
165 (*ostr) << "#-------------------------------------------------------------------------" << endl;
166 (*ostr) << "# You are running the CDF MidPoint plugin for FastJet " << endl;
167 (*ostr) << "# This is based on an implementation provided by Joey Huston. " << endl;
168 (*ostr) << "# If you use this plugin, please cite " << endl;
169 (*ostr) << "# G. C. Blazey et al., hep-ex/0005012. " << endl;
170 (*ostr) << "# in addition to the usual FastJet reference. " << endl;
171 (*ostr) << "#-------------------------------------------------------------------------" << endl;
172
173 // make sure we really have the output done.
174 ostr->flush();
175 }
176
177 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
178