fastjet 2.4.5
Classes | Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes
fastjet::DnnPlane Class Reference

class derived from DynamicNearestNeighbours that provides an implementation for the Euclidean plane More...

#include <DnnPlane.hh>

Inheritance diagram for fastjet::DnnPlane:
Inheritance graph
[legend]
Collaboration diagram for fastjet::DnnPlane:
Collaboration graph
[legend]

List of all members.

Classes

struct  SuperVertex
 Structure containing a vertex_handle and cached information on the nearest neighbour. More...

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 std::vector indices_to_remove, and add the points specified by the std::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
 CGAL object for dealing with triangulations.
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.
void _CrashIfVertexPresent (const Vertex_handle &vertex, const int &its_index)
 given a vertex_handle returned by CGAL on insertion of a new points, crash if it turns out that it corresponds to a vertex that we already knew about (usually because two points coincide)

Private Attributes

std::vector< SuperVertex_supervertex
bool _verbose
Triangulation _TR

Static Private Attributes

static const bool _crash_on_coincidence = true

Detailed Description

class derived from DynamicNearestNeighbours that provides an implementation for the Euclidean plane

Definition at line 45 of file DnnPlane.hh.


Constructor & Destructor Documentation

fastjet::DnnPlane::DnnPlane ( ) [inline]

empty initaliser

Definition at line 48 of file DnnPlane.hh.

{}
fastjet::DnnPlane::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.


Member Function Documentation

void fastjet::DnnPlane::_CrashIfVertexPresent ( const Vertex_handle vertex,
const int &  its_index 
) [private]

given a vertex_handle returned by CGAL on insertion of a new points, crash if it turns out that it corresponds to a vertex that we already knew about (usually because two points coincide)

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 fastjet::NEW_VERTEX.

                                                             {
  if (!_crash_on_coincidence) return;

  // vertices that do not have the same geometric position as any
  // other vertex so far added have info().val() == NEW_VERTEX -- this
  // is ensured by the InitialisedInt class, which forms the "info"
  // part of our
  // CGAL::Triangulation_vertex_base_with_info_2<InitialisedInt,K>
  //
  // If the vertex coincides with one that already exists, then
  // info().val() it's info().val() will have been updated (in
  // DNN:DNN) to be equal to a vertex "index".
  if (vertex->info().val() != NEW_VERTEX) {
    ostringstream err;
    err << "ERROR in DnnPlane::_CrashIfVertexPresent"
         <<endl << "Point "<<its_index<<" coincides with point "
         <<vertex->info().val() << endl;
    throw DnnError(err.str());
  } 
}
double fastjet::DnnPlane::_euclid_distance ( const Point p1,
const Point p2 
) const [inline, private]

CGAL object for dealing with triangulations.

calculates and returns the euclidean distance between points p1 and p2

Definition at line 102 of file DnnPlane.hh.

References fastjet::Point::x, and fastjet::Point::y.

                                                                         {
    double distx= p1.x()-p2.x();
    double disty= p1.y()-p2.y();
    return distx*distx+disty*disty;
  }
void fastjet::DnnPlane::_SetAndUpdateNearest ( const int &  j,
std::vector< int > &  indices_of_updated_neighbours 
) [private]

Determines and stores the nearest neighbour of j.

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.

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.

NB: note that we have _SetAndUpdateNearest as a completely separate routine from _SetNearest because we want to use one single ciruclation over voronoi neighbours to find the nearest neighbour and to update the voronoi neighbours if need be.

Definition at line 311 of file DnnPlane.cc.

References fastjet::HUGE_DOUBLE, and fastjet::INFINITE_VERTEX.

                                                                       {

  Vertex_handle current = _supervertex[j].vertex;
  Vertex_circulator vc = _TR.incident_vertices(current);
  Vertex_circulator done = vc;
  double dist;
  double mindist = HUGE_DOUBLE; // change this to "HUGE" or max_double?
  Vertex_handle nearest = _TR.infinite_vertex();

  // when there is only one finite point left in the triangulation, 
  // there are no triangles. Presumably this is why voronoi returns
  // NULL for the incident vertex circulator. Check if this is
  // happening before circulating over it... (Otherwise it crashes
  // when looking for neighbours of last point)
  if (vc != NULL) do { 
    if (vc->info().val() != INFINITE_VERTEX) {
      if (_verbose) cout << current->info().val() << " " << vc->info().val() << endl;
      // find distance between j and its Voronoi neighbour (vc)
      dist = _euclid_distance(current->point(), vc->point());
      // update the mindist if we are closer than anything found so far
      if (dist < mindist) {
        mindist = dist; nearest = vc; 
        if (_verbose) cout << "nearer ";
      } 
      // find index corresponding to vc for easy manipulation
      int vcindx = vc->info().val();
      if (_verbose) cout << vc->point() << "; "<< dist << endl;
      // check if j is closer to vc than vc's currently registered
      // nearest neighbour (and update things if it is)
      if (dist < _supervertex[vcindx].NNdistance) {
        if (_verbose) cout << vcindx << "'s NN becomes " << current->info().val() << endl;
        _supervertex[vcindx].NNdistance = dist;
        _supervertex[vcindx].NNindex = j;
        indices_of_updated_neighbours.push_back(vcindx);
      }
    }
  } while (++vc != done); // move on to next Voronoi neighbour
  // set j's supervertex info about nearest neighbour
  _supervertex[j].NNindex = nearest->info().val();
  _supervertex[j].NNdistance = mindist;
}
void fastjet::DnnPlane::_SetNearest ( const int &  j) [private]

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 fastjet::HUGE_DOUBLE, and fastjet::INFINITE_VERTEX.

                                         {
  Vertex_handle current = _supervertex[j].vertex;
  Vertex_circulator vc = _TR.incident_vertices(current);
  Vertex_circulator done = vc;
  double dist;
  double mindist = HUGE_DOUBLE; // change this to "HUGE" or max_double?
  Vertex_handle nearest = _TR.infinite_vertex();

  // when there is only one finite point left in the triangulation, 
  // there are no triangles. Presumably this is why voronoi returns
  // NULL for the incident vertex circulator. Check if this is
  // happening before circulating over it... (Otherwise it crashes
  // when looking for neighbours of last point)
  if (vc != NULL) do { 
    if ( vc->info().val() != INFINITE_VERTEX) {
      // find distance between j and its Voronoi neighbour (vc)
      if (_verbose) cout << current->info().val() << " " << vc->info().val() << endl;
      dist = _euclid_distance(current->point(), vc->point());
      // check if j is closer to vc than vc's currently registered
      // nearest neighbour (and update things if it is)
      if (dist < mindist) {
        mindist = dist; nearest = vc; 
        if (_verbose) cout << "nearer ";
      } 
      if (_verbose) cout << vc->point() << "; "<< dist << endl;
    }
  } while (++vc != done); // move on to next Voronoi neighbour
  // set j's supervertex info about nearest neighbour
  _supervertex[j].NNindex = nearest->info().val();
  _supervertex[j].NNdistance = mindist;
}
double fastjet::DnnPlane::eta ( const int  i) const [inline]

returns the eta point with index i.

Definition at line 152 of file DnnPlane.hh.

                                             {
  return _supervertex[i].vertex->point().x(); }
EtaPhi fastjet::DnnPlane::etaphi ( const int  i) const [inline]

returns the EtaPhi of point with index i.

Definition at line 148 of file DnnPlane.hh.

References fastjet::Point::x, and fastjet::Point::y.

                                                {
  Point * p = & (_supervertex[i].vertex->point());
  return EtaPhi(p->x(),p->y()); }
double fastjet::DnnPlane::NearestNeighbourDistance ( const int &  ii) const [inline, virtual]

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.

                                                                     {
  return _supervertex[ii].NNdistance;}
int fastjet::DnnPlane::NearestNeighbourIndex ( const int &  ii) const [inline, virtual]

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.

                                                               {
  return _supervertex[ii].NNindex;}
double fastjet::DnnPlane::phi ( const int  i) const [inline]

returns the phi point with index i.

Definition at line 155 of file DnnPlane.hh.

                                             {
  return _supervertex[i].vertex->point().y(); }
void fastjet::DnnPlane::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 
) [virtual]

remove the points labelled by the std::vector indices_to_remove, and add the points specified by the std::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.

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.

insertion and removal of points

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.

bool fastjet::DnnPlane::Valid ( const int &  index) const [inline, virtual]

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.

                                                   {
  if (index >= 0 && index < static_cast<int>(_supervertex.size())) {
    return (_supervertex[index].vertex != NULL);} else {return false;} }

Member Data Documentation

const bool fastjet::DnnPlane::_crash_on_coincidence = true [static, private]

Definition at line 95 of file DnnPlane.hh.

Definition at line 91 of file DnnPlane.hh.

Definition at line 98 of file DnnPlane.hh.

Definition at line 93 of file DnnPlane.hh.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines