#include <DnnPlane.hh>
Inheritance diagram for fastjet::DnnPlane:
Public Member Functions | |
DnnPlane () | |
empty initaliser | |
DnnPlane (const std::vector< EtaPhi > &, const bool &verbose=false) | |
Initialiser from a set of points on an Eta-Phi plane, where both eta and phi can have arbitrary ranges. | |
int | NearestNeighbourIndex (const int &ii) const |
Returns the index of the nearest neighbour of point labelled by ii (assumes ii is valid). | |
double | NearestNeighbourDistance (const int &ii) const |
Returns the distance to the nearest neighbour of point labelled by index ii (assumes ii is valid). | |
bool | Valid (const int &index) const |
Returns true iff the given index corresponds to a point that exists in the DNN structure (meaning that it has been added, and not removed in the meantime). | |
void | RemoveAndAddPoints (const std::vector< int > &indices_to_remove, const std::vector< EtaPhi > &points_to_add, std::vector< int > &indices_added, std::vector< int > &indices_of_updated_neighbours) |
remove the points labelled by the vector indices_to_remove, and add the points specified by the vector points_to_add (corresponding indices will be calculated automatically); the idea behind this routine is that the points to be added will somehow be close to the one or other of the points being removed and this can be used by the implementation to provide hints for inserting the new points in whatever structure it is using. | |
EtaPhi | etaphi (const int i) const |
returns the EtaPhi of point with index i. | |
double | eta (const int i) const |
returns the eta point with index i. | |
double | phi (const int i) const |
returns the phi point with index i. | |
Private Member Functions | |
double | _euclid_distance (const Point &p1, const Point &p2) const |
calculates and returns the euclidean distance between points p1 and p2 | |
void | _SetNearest (const int &j) |
Determines the index and distance of the nearest neighbour to point j and puts the information into the _supervertex entry for j. | |
void | _SetAndUpdateNearest (const int &j, std::vector< int > &indices_of_updated_neighbours) |
Determines and stores the nearest neighbour of j, and where necessary update the nearest-neighbour info of Voronoi neighbours of j;. | |
void | _CrashIfVertexPresent (const Vertex_handle &vertex, const int &its_index) |
Crashes if the given vertex handle already exists. | |
Private Attributes | |
std::vector< SuperVertex > | _supervertex |
bool | _verbose |
Triangulation | _TR |
Static Private Attributes | |
static const bool | _crash_on_coincidence = true |
Classes | |
struct | SuperVertex |
Structure containing a vertex_handle and cached information on the nearest neighbour. More... |
Definition at line 45 of file DnnPlane.hh.
|
empty initaliser
Definition at line 48 of file DnnPlane.hh. 00048 {}
|
|
Initialiser from a set of points on an Eta-Phi plane, where both eta and phi can have arbitrary ranges.
Definition at line 44 of file DnnPlane.cc. References _CrashIfVertexPresent(), _SetNearest(), _supervertex, _TR, _verbose, fastjet::INFINITE_VERTEX, and fastjet::DnnPlane::SuperVertex::vertex. 00045 { 00046 00047 _verbose = verbose; 00048 int n = input_points.size(); 00049 00050 // construct Voronoi diagram in such a way as to get the vertex handles 00051 // and remember to set CGAL info with the index of the vertex 00052 SuperVertex sv; 00053 for (int i = 0; i < n; i++) { 00054 sv.vertex = 00055 _TR.insert(Point(input_points[i].first, input_points[i].second)); 00056 00057 // we are not up to dealing with coincident vertices, so make 00058 // sure the user knows! 00059 _CrashIfVertexPresent(sv.vertex, i); 00060 00061 // we need to assicate an index to each vertex -- thus when we get 00062 // a vertex (e.g. as a nearest neighbour) from CGAL, we will be 00063 // able to figure out which particle it corresponded to. 00064 sv.vertex->info() = i; 00065 _supervertex.push_back(sv); 00066 } 00067 00068 // label infinite vertex info with negative index 00069 _TR.infinite_vertex()->info() = INFINITE_VERTEX; 00070 00071 // set up the structure that holds nearest distances and neighbours 00072 for (int j = 0; j < n; j++) {_SetNearest(j);} 00073 00074 }
|
|
Crashes if the given vertex handle already exists. Otherwise it does the bookkeeping for future such tests Definition at line 80 of file DnnPlane.cc. References _crash_on_coincidence, and fastjet::NEW_VERTEX. Referenced by DnnPlane(), and RemoveAndAddPoints(). 00081 { 00082 if (!_crash_on_coincidence) return; 00083 00084 // vertices that do not have the same geometric position as any 00085 // other vertex so far added have info().val() == NEW_VERTEX -- this 00086 // is ensured by the InitialisedInt class, which forms the "info" 00087 // part of our 00088 // CGAL::Triangulation_vertex_base_with_info_2<InitialisedInt,K> 00089 // 00090 // If the vertex coincides with one that already exists, then 00091 // info().val() it's info().val() will have been updated (in 00092 // DNN:DNN) to be equal to a vertex "index". 00093 if (vertex->info().val() != NEW_VERTEX) { 00094 ostringstream err; 00095 err << "ERROR in DnnPlane::_CrashIfVertexPresent" 00096 <<endl << "Point "<<its_index<<" coincides with point " 00097 <<vertex->info().val() << endl; 00098 throw DnnError(err.str()); 00099 } 00100 }
|
|
calculates and returns the euclidean distance between points p1 and p2
Definition at line 102 of file DnnPlane.hh. Referenced by _SetAndUpdateNearest(), and _SetNearest(). 00102 { 00103 double distx= p1.x()-p2.x(); 00104 double disty= p1.y()-p2.y(); 00105 return distx*distx+disty*disty; 00106 }
|
|
Determines and stores the nearest neighbour of j, and where necessary update the nearest-neighbour info of Voronoi neighbours of j;. For each voronoi neighbour D of j if the distance between j and D is less than D's own nearest neighbour, then update the nearest-neighbour info in D; push D's index onto indices_of_updated_neighbours Note that j is NOT pushed onto indices_of_updated_neighbours -- if you want it there, put it there yourself. Definition at line 311 of file DnnPlane.cc. References _euclid_distance(), _supervertex, _TR, _verbose, fastjet::HUGE_DOUBLE, and fastjet::INFINITE_VERTEX. Referenced by RemoveAndAddPoints(). 00313 { 00314 00315 Vertex_handle current = _supervertex[j].vertex; 00316 Vertex_circulator vc = _TR.incident_vertices(current); 00317 Vertex_circulator done = vc; 00318 double dist; 00319 double mindist = HUGE_DOUBLE; // change this to "HUGE" or max_double? 00320 Vertex_handle nearest = _TR.infinite_vertex(); 00321 00322 // when there is only one finite point left in the triangulation, 00323 // there are no triangles. Presumably this is why voronoi returns 00324 // NULL for the incident vertex circulator. Check if this is 00325 // happening before circulating over it... (Otherwise it crashes 00326 // when looking for neighbours of last point) 00327 if (vc != NULL) do { 00328 if (vc->info().val() != INFINITE_VERTEX) { 00329 if (_verbose) cout << current->info().val() << " " << vc->info().val() << endl; 00330 // find distance between j and its Voronoi neighbour (vc) 00331 dist = _euclid_distance(current->point(), vc->point()); 00332 // update the mindist if we are closer than anything found so far 00333 if (dist < mindist) { 00334 mindist = dist; nearest = vc; 00335 if (_verbose) cout << "nearer "; 00336 } 00337 // find index corresponding to vc for easy manipulation 00338 int vcindx = vc->info().val(); 00339 if (_verbose) cout << vc->point() << "; "<< dist << endl; 00340 // check if j is closer to vc than vc's currently registered 00341 // nearest neighbour (and update things if it is) 00342 if (dist < _supervertex[vcindx].NNdistance) { 00343 if (_verbose) cout << vcindx << "'s NN becomes " << current->info().val() << endl; 00344 _supervertex[vcindx].NNdistance = dist; 00345 _supervertex[vcindx].NNindex = j; 00346 indices_of_updated_neighbours.push_back(vcindx); 00347 } 00348 } 00349 } while (++vc != done); // move on to next Voronoi neighbour 00350 // set j's supervertex info about nearest neighbour 00351 _supervertex[j].NNindex = nearest->info().val(); 00352 _supervertex[j].NNdistance = mindist; 00353 }
|
|
Determines the index and distance of the nearest neighbour to point j and puts the information into the _supervertex entry for j.
Definition at line 261 of file DnnPlane.cc. References _euclid_distance(), _supervertex, _TR, _verbose, fastjet::HUGE_DOUBLE, and fastjet::INFINITE_VERTEX. Referenced by DnnPlane(), and RemoveAndAddPoints(). 00261 { 00262 Vertex_handle current = _supervertex[j].vertex; 00263 Vertex_circulator vc = _TR.incident_vertices(current); 00264 Vertex_circulator done = vc; 00265 double dist; 00266 double mindist = HUGE_DOUBLE; // change this to "HUGE" or max_double? 00267 Vertex_handle nearest = _TR.infinite_vertex(); 00268 00269 // when there is only one finite point left in the triangulation, 00270 // there are no triangles. Presumably this is why voronoi returns 00271 // NULL for the incident vertex circulator. Check if this is 00272 // happening before circulating over it... (Otherwise it crashes 00273 // when looking for neighbours of last point) 00274 if (vc != NULL) do { 00275 if ( vc->info().val() != INFINITE_VERTEX) { 00276 // find distance between j and its Voronoi neighbour (vc) 00277 if (_verbose) cout << current->info().val() << " " << vc->info().val() << endl; 00278 dist = _euclid_distance(current->point(), vc->point()); 00279 // check if j is closer to vc than vc's currently registered 00280 // nearest neighbour (and update things if it is) 00281 if (dist < mindist) { 00282 mindist = dist; nearest = vc; 00283 if (_verbose) cout << "nearer "; 00284 } 00285 if (_verbose) cout << vc->point() << "; "<< dist << endl; 00286 } 00287 } while (++vc != done); // move on to next Voronoi neighbour 00288 // set j's supervertex info about nearest neighbour 00289 _supervertex[j].NNindex = nearest->info().val(); 00290 _supervertex[j].NNdistance = mindist; 00291 }
|
|
returns the eta point with index i.
Definition at line 152 of file DnnPlane.hh. References _supervertex. 00152 { 00153 return _supervertex[i].vertex->point().x(); }
|
|
returns the EtaPhi of point with index i.
Definition at line 148 of file DnnPlane.hh. References _supervertex. Referenced by fastjet::Dnn2piCylinder::_CreateNecessaryMirrorPoints(). 00148 { 00149 Point * p = & (_supervertex[i].vertex->point()); 00150 return EtaPhi(p->x(),p->y()); }
|
|
Returns the distance to the nearest neighbour of point labelled by index ii (assumes ii is valid).
Implements fastjet::DynamicNearestNeighbours. Definition at line 141 of file DnnPlane.hh. References _supervertex. Referenced by fastjet::Dnn2piCylinder::_CreateNecessaryMirrorPoints(), fastjet::Dnn4piCylinder::NearestNeighbourDistance(), fastjet::Dnn3piCylinder::NearestNeighbourDistance(), fastjet::Dnn2piCylinder::NearestNeighbourDistance(), fastjet::Dnn4piCylinder::NearestNeighbourIndex(), fastjet::Dnn3piCylinder::NearestNeighbourIndex(), and fastjet::Dnn2piCylinder::NearestNeighbourIndex(). 00141 { 00142 return _supervertex[ii].NNdistance;}
|
|
Returns the index of the nearest neighbour of point labelled by ii (assumes ii is valid).
Implements fastjet::DynamicNearestNeighbours. Definition at line 138 of file DnnPlane.hh. References _supervertex. Referenced by fastjet::Dnn4piCylinder::NearestNeighbourIndex(), fastjet::Dnn3piCylinder::NearestNeighbourIndex(), and fastjet::Dnn2piCylinder::NearestNeighbourIndex(). 00138 { 00139 return _supervertex[ii].NNindex;}
|
|
returns the phi point with index i.
Definition at line 155 of file DnnPlane.hh. References _supervertex. 00155 { 00156 return _supervertex[i].vertex->point().y(); }
|
|
remove the points labelled by the vector indices_to_remove, and add the points specified by the vector points_to_add (corresponding indices will be calculated automatically); the idea behind this routine is that the points to be added will somehow be close to the one or other of the points being removed and this can be used by the implementation to provide hints for inserting the new points in whatever structure it is using. In a kt-algorithm the points being added will be a result of a combination of the points to be removed -- hence the proximity is (more or less) guaranteed. Implements fastjet::DynamicNearestNeighbours. Definition at line 114 of file DnnPlane.cc. References _CrashIfVertexPresent(), _SetAndUpdateNearest(), _SetNearest(), _supervertex, _TR, _verbose, and fastjet::INFINITE_VERTEX. Referenced by fastjet::Dnn2piCylinder::_CreateNecessaryMirrorPoints(), fastjet::Dnn4piCylinder::RemoveAndAddPoints(), fastjet::Dnn3piCylinder::RemoveAndAddPoints(), and fastjet::Dnn2piCylinder::RemoveAndAddPoints(). 00118 { 00119 00120 00121 // build set of UNION of Voronoi neighbours of a pair of nearest 00122 // neighbours 00123 set<int> NeighbourUnion; 00124 // later on it will be convenient to have access to a set (rather 00125 // than vector) of indices being removed 00126 set<int> indices_removed; 00127 00128 // for each of the indices to be removed add the voronoi neighbourhood to 00129 // the NeighbourUnion set. 00130 for (size_t ir = 0; ir < indices_to_remove.size(); ir++) { 00131 int index = indices_to_remove[ir]; 00132 indices_removed.insert(index); 00133 if (_verbose) cout << " Starting RemoveAndAddPoints" << endl; 00134 if (_verbose) cout << " point " << index << endl; 00135 // have a circulators that will go round the Voronoi neighbours of 00136 // _supervertex[index1].vertex 00137 Vertex_circulator vc = _TR.incident_vertices(_supervertex[index].vertex); 00138 Vertex_circulator done = vc; 00139 do { 00140 // if a neighbouring vertex not the infinite vertex, then add it 00141 // to our union of neighbouring vertices. 00142 if (_verbose) cout << "examining " << vc->info().val() << endl; 00143 if (vc->info().val() != INFINITE_VERTEX) { 00144 // NB: from it=1 onwards occasionally it might already have 00145 // been inserted -- but double insertion still leaves only one 00146 // copy in the set, so there's no problem 00147 NeighbourUnion.insert(vc->info().val()); 00148 if (_verbose) cout << "inserted " << vc->info().val() << endl; 00149 } 00150 } while (++vc != done); 00151 } 00152 00153 if (_verbose) { 00154 set<int>::iterator it = NeighbourUnion.begin(); 00155 cout << "Union of neighbours of combined points" << endl; 00156 for ( ; it != NeighbourUnion.end(); ++it ) { 00157 cout << *it << endl; 00158 } 00159 } 00160 00161 // update set, triangulation and supervertex info 00162 for (size_t ir = 0; ir < indices_to_remove.size(); ir++) { 00163 int index = indices_to_remove[ir]; 00164 00165 // NeighbourUnion should not contain the points to be removed 00166 // (because later we will assume they still exist). 00167 NeighbourUnion.erase(indices_to_remove[ir]); 00168 00169 // points to be removed should also be eliminated from the 00170 // triangulation and the supervertex structure should be updated 00171 // to reflect the fact that the points are no longer valid. 00172 _TR.remove(_supervertex[index].vertex); 00173 _supervertex[index].vertex = NULL; 00174 } 00175 00176 // add new point: give a "hint" to the inserter that 00177 // the new point should be added close to old points -- the easiest way 00178 // of getting this is to take a point from the NeighbourUnion AFTER we have 00179 // removed point1, point2, and to get one of its incident faces. 00180 // 00181 // This hinting improves speed by c. 25% for 10^4 points because it 00182 // avoids the need for a costly (sqrt{N}) location step (at least 00183 // with a non-hierarchical triangulation -- with a hierarchical one, 00184 // this step could be done away with, though there will still be a 00185 // cost of O(ln N) to pay. 00186 // 00187 // For some reason inserting the point before the two removals 00188 // slows things down by c. 25%. This importance of the order 00189 // is not understood. 00190 // 00191 // At some point it might be worth trying to select the "nearest" 00192 // of the various points in the neighbour union to avoid large 00193 // steps in cases where we have 0..2pi periodicity and the first member 00194 // of the neighbour union happens to be on the wrong side. 00195 Face_handle face; 00196 if (indices_to_remove.size() > 0) { 00197 // face can only be found if there were points to remove in first place 00198 face = _TR.incident_faces( 00199 _supervertex[*NeighbourUnion.begin()].vertex);} 00200 // make sure the output arrays are empty 00201 indices_added.clear(); 00202 indices_of_updated_neighbours.clear(); 00203 for (size_t ia = 0; ia < points_to_add.size(); ia++) { 00204 SuperVertex sv; 00205 _supervertex.push_back(sv); 00206 int index = _supervertex.size()-1; 00207 indices_added.push_back(index); 00208 00209 if (indices_to_remove.size() > 0) { 00210 // be careful of using face (for location hinting) only when it exists 00211 _supervertex[index].vertex = _TR.insert(Point(points_to_add[ia].first, 00212 points_to_add[ia].second),face);} 00213 else { 00214 _supervertex[index].vertex = _TR.insert(Point(points_to_add[ia].first, 00215 points_to_add[ia].second)); 00216 } 00217 // we are not up to dealing with coincident vertices, so make 00218 // sure the user knows! 00219 _CrashIfVertexPresent(_supervertex[index].vertex, index); 00220 _supervertex[index].vertex->info() = index; 00221 00222 // first find nearest neighbour of "newpoint" (shorthand for 00223 // _supervertex[index].vertex); while we're at it, for each of the 00224 // voronoi neighbours, "D", of newpoint, examine whether newpoint is 00225 // closer to "D" than D's current nearest neighbour -- when this 00226 // occurs, put D into indices_of_updated_neighbours. 00227 // 00228 // manually put newpoint on indices_of_updated_neighbours 00229 indices_of_updated_neighbours.push_back(index); 00230 _SetAndUpdateNearest(index, indices_of_updated_neighbours); 00231 } 00232 00233 // for Voronoi neighbours j of any of the removed points for which 00234 // one of those removed points was the nearest neighbour, 00235 // redetermine the nearest neighbour of j and add j onto the vector 00236 // of indices_of_updated_neighbours. 00237 set<int>::iterator it2 = NeighbourUnion.begin(); 00238 for ( ; it2 != NeighbourUnion.end(); ++it2 ) { 00239 int j = *it2; 00240 // the if avoids the vertex at infinity, which gets a negative index 00241 if( j != INFINITE_VERTEX ) { 00242 // this is where we check if the nearest neighbour of j was one 00243 // of the removed points 00244 if (indices_removed.count(_supervertex[j].NNindex)) { 00245 if (_verbose) cout << "j " << j << endl; 00246 _SetNearest(j); 00247 indices_of_updated_neighbours.push_back(j); 00248 if (_verbose) cout << "NN of " << j << " : " 00249 << _supervertex[j].NNindex 00250 << ", dist = " << _supervertex[j].NNdistance <<endl; 00251 } 00252 } 00253 } 00254 00255 }
|
|
Returns true iff the given index corresponds to a point that exists in the DNN structure (meaning that it has been added, and not removed in the meantime).
Implements fastjet::DynamicNearestNeighbours. Definition at line 144 of file DnnPlane.hh. References _supervertex. Referenced by fastjet::Dnn4piCylinder::Valid(), fastjet::Dnn3piCylinder::Valid(), and fastjet::Dnn2piCylinder::Valid(). 00144 { 00145 if (index >= 0 && index < static_cast<int>(_supervertex.size())) { 00146 return (_supervertex[index].vertex != NULL);} else {return false;} }
|
|
Definition at line 95 of file DnnPlane.hh. Referenced by _CrashIfVertexPresent(). |
|
Definition at line 91 of file DnnPlane.hh. Referenced by _SetAndUpdateNearest(), _SetNearest(), DnnPlane(), eta(), etaphi(), NearestNeighbourDistance(), NearestNeighbourIndex(), phi(), RemoveAndAddPoints(), and Valid(). |
|
Definition at line 98 of file DnnPlane.hh. Referenced by _SetAndUpdateNearest(), _SetNearest(), DnnPlane(), and RemoveAndAddPoints(). |
|
Definition at line 93 of file DnnPlane.hh. Referenced by _SetAndUpdateNearest(), _SetNearest(), DnnPlane(), and RemoveAndAddPoints(). |