Skip to Content.
Sympa Menu

cgal-discuss - Re: [cgal-discuss] How to obtain the index of the point in dD spatial search

Subject: CGAL users discussion list

List archive

Re: [cgal-discuss] How to obtain the index of the point in dD spatial search


Chronological Thread 
  • From: "Al Hasan" <>
  • To:
  • Subject: Re: [cgal-discuss] How to obtain the index of the point in dD spatial search
  • Date: Thu, 15 May 2008 14:55:06 -0400
  • Domainkey-signature: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=message-id:date:from:to:subject:in-reply-to:mime-version:content-type:references; b=HTpMCOYQsVwQ6K57SztfgMUdgiPIgA+45m8yR+be1PDt75FwVN1xNvn/BkNbaUPmKmroRM2/FtsA1xV1gs1N88o1Xc7ywQeBVbP/9x9cj6oZZwbm7vNfln9poIgrt2bKdspXT11LwswS3uy3yZzgmttqjTYlk6CqZg1zDTekdpw=

You can write your own customized point and distance class. In the poing
class you can make a member variable (index of the point), . It will provide
an
way to do what you want. Read "user defined point and distance" section
in CGAL manual.

I have done it and it is working fine. I am attaching my code for that.

Hasan

On Sun, May 11, 2008 at 7:04 PM, Jian Sun
<>
wrote:
> Hi Kevin,
>
> Thank you for your suggestion. It is very helpful. But it does not works as
> it is.
> One needs to define a comparison operator on Point.
>
> It serves my purpose well now, though it is different from the solution that
> I am looking for.
> I am looking for a way to customize Point_d as one does for the vertex in
> the Delaunay
> triangulation. If that is doable, we can add an index into each point, which
> avoids the
> extra step of looking up in map.
>
> sunjian
>
> On Sun, May 11, 2008 at 2:12 AM, Kevin Xu
> <>
> wrote:
>>
>> Hi, Sun Jian,
>>
>> You can use the map container provided by STL to achieve this.
>>
>> Let's take the range search (with kd-tree) for example:
>>
>> KD_tree STree; // define a kd tree for all vertices
>> std::map<Point,int> mVPToID; // map vertex pos. to index
>> std::map<Point,int>::iterator iPTI; // an iterator for mVPToID
>>
>> // construct mapping
>> for (all vertices in a mesh or other data structure)
>> {
>> STree.insert (VP(i)); // add position of the ith vertex
>> mVPToID.insert (std::make_pair(VP(i), ID(i))); // make a map
>> between pos. & id
>> }
>>
>> // perform range searching
>> std::vector<Point> FP; // store all found points (pos.)
>> Fuzzy_sphere fs(query_pos, dRadia, 0.0); // a fuzzy query sphere
>> STree.search(std::back_inserter(FP), fs); // do searching
>>
>> // collect ids for the found points
>> for (int i=0; i<FP.size(); i++)
>> {
>> if ((iPTI=mVPToID.find(FP[i])) != mVPToID.end())
>> { // find the record of FP[i] to get its id
>> ID(i) = iPTI->second; // id of the ith point is "iPTI->second"
>> }
>> }
>>
>> /// end of the routine
>>
>> The routines for k-NN search, etc., are similar.
>>
>> HTH.
>>
>> Good luck!
>> --
>> You are currently subscribed to cgal-discuss.
>> To unsubscribe or access the archives, go to
>> https://lists-sop.inria.fr/wws/info/cgal-discuss
>
>



--
--------------------------------
Mohammad Al Hasan
2257 Burdett Ave,
Troy NY 12180
template <unsigned int D>
struct Distance {
  typedef Point<D> Point_D;
  typedef Point<D> Query_item;

  double transformed_distance(const Point_D& p1, const Point_D& p2) const {
    double ret_val(0.0);
    for (int i=0; i<D; i++) {
      double tempval = p1.get_index_val(i) - p2.get_index_val(i);
      ret_val += tempval*tempval;
    }
    return ret_val;
  }

  template <class TreeTraits>
  double min_distance_to_rectangle(const Point_D& p,
				   const CGAL::Kd_tree_rectangle<TreeTraits>& b) const {   
    double distance(0.0); 
    for (int i=0; i<D; i++) {
      double h = p.get_index_val(i);
      if (h < b.min_coord(i)) distance += (b.min_coord(i)-h)*(b.min_coord(i)-h);
      if (h > b.max_coord(i)) distance += (h-b.max_coord(i))*(h-b.max_coord(i));
    }
    return distance;
  }
  
  template <class TreeTraits>
  double max_distance_to_rectangle(const Point_D& p,
				   const CGAL::Kd_tree_rectangle<TreeTraits>& b) const {   

    double distance(0.0);
    for (int i=0; i<D; i++) {
      double h = p.get_index_val(i);
      double tempval = (h >= (b.min_coord(i)+b.max_coord(i))/2.0) ?
                (h-b.min_coord(i))*(h-b.min_coord(i)) : (b.max_coord(i)-h)*(b.max_coord(i)-h);
      distance += tempval;
    }
    return distance;
    
  }
  
  double new_distance(double& dist, double old_off, double new_off,
		      int cutting_dimension)  const {
    return dist + new_off*new_off - old_off*old_off;
  }
  
  double transformed_distance(double d) const { return d*d; }
  
  double inverse_of_transformed_distance(double d) { return std::sqrt(d); }
  
}; // end of struct Distance
#include <iostream>

using namespace std;

typedef unsigned int UINT;

template<UINT MAXSIZE> struct Point;

template<UINT MAXSIZE>
ostream& operator<< (ostream&, const Point<MAXSIZE>& p);


template <UINT MAXSIZE>
struct Point {
  int id;
  double vec[MAXSIZE];
  Point() { 
    for (int i=0; i<MAXSIZE; i++)
      vec[i]= 0;
  }
  unsigned int dimension() const { return MAXSIZE;}

  const double * get_data() const { return vec;}
  Point(int id, double input[]) {
    this->id = id;
    for (int i=0; i<MAXSIZE; i++)
      vec[i]= input[i];
  }

  double get_index_val(int i)const {return vec[i];}

  friend ostream& operator << <>(ostream&, const Point<MAXSIZE>&); 

  bool operator==(const Point& p) const 
  {
    for (int i=0; i<MAXSIZE; i++)
      if (vec[i] != p.vec[i]) return false;
    return true;
  }

  bool  operator!=(const Point& p) const { return ! (*this == p); }

};
template<unsigned int MAXSIZE>
ostream& operator<< (ostream& ostr, const Point<MAXSIZE>& p) {
  ostr << p.id;
}
//file: examples/Spatial_searching/User_defined_point_and_distance.C

#include <CGAL/basic.h>
#include <CGAL/Search_traits.h>
#include <CGAL/Orthogonal_k_neighbor_search.h>
#include <CGAL/Fuzzy_sphere.h>
#include "MyPoint.h"  // defines types Point, Construct_coord_iterator
#include "MyDistance.h"

typedef Point<5> Point_D;
typedef Distance<5> Distance_D;

namespace CGAL {

  template <>
  struct Kernel_traits<Point_D> {
    struct Kernel {
      typedef double FT;
      typedef double RT;
    };
  };
}

struct Construct_coord_iterator {
  const double* operator()(const Point_D& p) const 
  { return static_cast<const double*>(p.vec); }

  const double* operator()(const Point_D& p, int)  const
  { 
    int dim = p.dimension(); 
    return static_cast<const double*>(p.vec+dim); 
  }
};


typedef CGAL::Search_traits<double, Point_D, const double*, Construct_coord_iterator> Traits;
typedef CGAL::Fuzzy_sphere<Traits> Fuzzy_sphere;
typedef CGAL::Orthogonal_k_neighbor_search<Traits, Distance_D > K_neighbor_search;
typedef K_neighbor_search::Tree Tree;

int main() {
  const int N = 1000;
  const int K = 3;
  // generator for random data points in the cube ( (-1,-1,-1), (1,1,1) ) 
  // Random_points_iterator rpit( 5);
  std::list<Point_D> points; 
  double p1d[5] = {1,0,0,0,0};Point_D p1(1, p1d);
  double p2d[5] = {0,2,0,0,0};Point_D p2(2, p2d);
  double p3d[5] = {0,0,3,0,0};Point_D p3(3, p3d);
  double p4d[5] = {0,0,0,4,0};Point_D p4(4, p4d);
  double p5d[5] = {0,0,0,0,5};Point_D p5(5, p5d);
  points.push_back(p1); points.push_back(p2); points.push_back(p3); points.push_back(p4); points.push_back(p5);
  // Insert number_of_data_points in the tree
  Tree tree(points.begin(),points.end());

  double query_data[5] = {1,0,0,0,0};
  Point_D query(0, query_data);
  Distance_D tr_dist;

  // search K nearest neighbours
  K_neighbor_search search(tree, query, K);
  for(K_neighbor_search::iterator it = search.begin(); it != search.end(); it++){
    std::cout << " d(q, nearest neighbor)=  " << it->first << ","
	      << tr_dist.inverse_of_transformed_distance(it->second) << std::endl; 
  }
  // search K furthest neighbour searching, with eps=0, search_nearest=false 
  K_neighbor_search search2(tree, query, K, 0.0, false);
  
  for(K_neighbor_search::iterator it = search2.begin(); it != search2.end(); it++){
    std::cout << " d(q, furthest neighbor)=  " 
	      << tr_dist.inverse_of_transformed_distance(it->second) << std::endl; 
  }
  Fuzzy_sphere exact_range(query, 3);
  std::list<Point_D> result;
  tree.search(std::back_inserter( result ), exact_range);
  std::cout << "The points in the circle centered at (0,0) is :";
  std::copy(result.begin(), result.end(), std::ostream_iterator<Point_D>(std::cout, " "));
  std::cout << std::endl;
  return 0;
}



Archive powered by MHonArc 2.6.16.

Top of Page