30 #include "fastjet/Error.hh"
31 #include "fastjet/PseudoJet.hh"
32 #include "fastjet/ClusterSequence.hh"
40 #ifndef DROP_CGAL // in case we do not have the code for CGAL
41 #include "fastjet/internal/Dnn4piCylinder.hh"
42 #include "fastjet/internal/Dnn3piCylinder.hh"
43 #include "fastjet/internal/Dnn2piCylinder.hh"
46 FASTJET_BEGIN_NAMESPACE
57 void ClusterSequence::_delaunay_cluster () {
61 vector<EtaPhi> points(n);
62 for (
int i = 0; i < n; i++) {
63 points[i] = EtaPhi(_jets[i].rap(),_jets[i].phi_02pi());
68 auto_ptr<DynamicNearestNeighbours> DNN;
69 #ifndef DROP_CGAL // strategy = NlnN* are not supported if we drop CGAL...
71 bool ignore_nearest_is_mirror = (_Rparam < twopi);
72 if (_strategy == NlnN4pi) {
73 DNN.reset(
new Dnn4piCylinder(points,verbose));
74 }
else if (_strategy == NlnN3pi) {
75 DNN.reset(
new Dnn3piCylinder(points,ignore_nearest_is_mirror,verbose));
76 }
else if (_strategy == NlnN) {
77 DNN.reset(
new Dnn2piCylinder(points,ignore_nearest_is_mirror,verbose));
80 if (_strategy == NlnN4pi || _strategy == NlnN3pi || _strategy == NlnN) {
82 err <<
"ERROR: Requested strategy "<<strategy_string()<<
" but it is not"<<endl;
83 err <<
" supported because FastJet was compiled without CGAL"<<endl;
84 throw Error(err.str());
90 err <<
"ERROR: Unrecognized value for strategy: "<<_strategy<<endl;
92 throw Error(err.str());
101 for (
int ii = 0; ii < n; ii++) {
102 _add_ktdistance_to_map(ii, DijMap, DNN.get());
108 for (
int i=0;i<n;i++) {
111 TwoVertices SmallestDijPair;
115 bool recombine_with_beam;
117 SmallestDij = DijMap.begin()->first;
118 SmallestDijPair = DijMap.begin()->second;
119 jet_i = SmallestDijPair.first;
120 jet_j = SmallestDijPair.second;
126 DijMap.erase(DijMap.begin());
131 recombine_with_beam = (jet_j == BeamJet);
132 if (!recombine_with_beam) {Valid2 = DNN->Valid(jet_j);}
133 else {Valid2 =
true;}
135 }
while ( !DNN->Valid(jet_i) || !Valid2);
141 if (! recombine_with_beam) {
143 _do_ij_recombination_step(jet_i, jet_j, SmallestDij, nn);
159 EtaPhi newpoint(_jets[nn].rap(), _jets[nn].phi_02pi());
161 points.push_back(newpoint);
164 _do_iB_recombination_step(jet_i, SmallestDij);
171 if (i == n-1) {
break;}
173 vector<int> updated_neighbours;
174 if (! recombine_with_beam) {
177 DNN->RemoveCombinedAddCombination(jet_i, jet_j,
178 points[points.size()-1], point3,
183 if (static_cast<unsigned int> (point3) != points.size()-1) {
184 throw Error(
"INTERNAL ERROR: point3 != points.size()-1");}
187 DNN->RemovePoint(jet_i, updated_neighbours);
191 vector<int>::iterator it = updated_neighbours.begin();
192 for (; it != updated_neighbours.end(); ++it) {
194 _add_ktdistance_to_map(ii, DijMap, DNN.get());
217 void ClusterSequence::_add_ktdistance_to_map(
220 const DynamicNearestNeighbours * DNN) {
222 double yiB = jet_scale_for_algorithm(_jets[ii]);
226 DijMap.insert(DijEntry(yiB, TwoVertices(ii,-1)));
228 double DeltaR2 = DNN->NearestNeighbourDistance(ii) * _invR2;
239 DijMap.insert(DijEntry(yiB, TwoVertices(ii,-1)));
241 double kt2i = jet_scale_for_algorithm(_jets[ii]);
242 int jj = DNN->NearestNeighbourIndex(ii);
243 if (kt2i <= jet_scale_for_algorithm(_jets[jj])) {
244 double dij = DeltaR2 * kt2i;
245 DijMap.insert(DijEntry(dij, TwoVertices(ii,jj)));
252 FASTJET_END_NAMESPACE