C++ Triangulation::locate方法代码示例

void test_case(int n) {

	// cout << "---- " << n << " ----" << endl;

	// Read all infected people
	std::vector<K::Point_2> infected;
	infected.reserve(n);

	for(int i=0; i<n; i++) {
		// cin >> x[i] >> y[i];
		K::Point_2 p;
		cin >> p;
		infected.push_back(p);
		cout << "Read in point " << p << endl;
	}

	// Construct Delauney triangulation
	Triangulation t;
	t.insert(infected.begin(), infected.end());

	// Read all healthy people
	int m;
	cin >> m;

	for(int i=0; i<m; i++) {
		K::Point_2 escaper;
		long d;
		cin >> escaper >> d;

		// --- Find an escape path for this person ---
		// Find out at which face we are
		Face_handle current_face = t.locate(escaper);

		// Check if we are already outside
		if(t.is_infinite(current_face)) {
			cout << "y";
			continue;
		}
		// Check if we are already getting infected
		/*K::Point_2 nearest_infected = t.nearest_vertex(escaper, current_face)->point();
		cout << "Nearest infected person: " << nearest_infected << endl;
		int dx = nearest_infected.x() - escaper.x();
		int dy = nearest_infected.y() - escaper.y();
		long nearest_sqd = dx * dx + dy * dy;
		if(nearest_sqd < d) {
			cout << "n";
			continue;
		}*/

		// Recurse
		vector<Face_handle> visited;
		bool result = recurse(current_face, d, visited, t);

		if(result)
			cout << "POSSIBLE TO ESCAPE" << endl;
		else
			cout << "COULDN'T ESCAPE :(" << endl;
	}

}

Face_handle test_point_location(const Triangulation &t,
                                const Point &query,
                                const Triangulation::Locate_type &lt_in)
{
  Triangulation::Locate_type lt, lt2;
  int li, li2;
  Face_handle fh;
  CGAL::Bounded_side bs;
  CGAL::Oriented_side os;

  fh = t.locate(query, lt, li);
  CGAL_assertion(lt == lt_in);
  if (lt_in == Triangulation::EMPTY) {
    CGAL_assertion(fh == Face_handle());
    return fh;
  }

  bs = t.side_of_face(query, fh, lt2, li2);
  os = t.oriented_side(fh, query);
  CGAL_USE(bs);
  CGAL_USE(os);
  
  CGAL_assertion(lt2 == lt_in);

  switch (lt_in)
    {
    case Triangulation::VERTEX:
    case Triangulation::EDGE:
    {
      CGAL_assertion(fh != Face_handle());
      CGAL_assertion(bs == CGAL::ON_BOUNDARY);
      CGAL_assertion(os == CGAL::ON_ORIENTED_BOUNDARY);

      CGAL_assertion(li == li2);
      break;
    }
    case Triangulation::FACE:
    {
      CGAL_assertion(fh != Face_handle());
      CGAL_assertion(bs == CGAL::ON_BOUNDED_SIDE);
      CGAL_assertion(os == CGAL::ON_POSITIVE_SIDE);
      break;
    }
    case Triangulation::EMPTY:
    {
      // Handled above
      CGAL_assertion(false);
      break;
    }
    case Triangulation::OUTSIDE_CONVEX_HULL: CGAL_error();
    case Triangulation::OUTSIDE_AFFINE_HULL: CGAL_error();
    }

  return fh;
}

bool
is_p_outside(const Point& q, const Triangulation& triang,
             const Cell_handle& start_cell,
             Cell_handle& c, int& u, int& v)
{
   Triangulation::Locate_type lt;
   u = -1; v = -1;
   c = triang.locate(q, lt, u, v, start_cell);
   if( lt == Triangulation::OUTSIDE_AFFINE_HULL )
   {
     cerr << "Point " << q << " is outside the affine hull." << endl;
     return true;
   }
   else if( lt == Triangulation::OUTSIDE_CONVEX_HULL )
     return true;
   else
   {
     if( lt == Triangulation::CELL )
     {
        if( c->outside )
           return true; 
        else
           return false; 
     }
     else if( lt == Triangulation::FACET )
     {
        Cell_handle _c = c->neighbor(u);
        if( c->outside && _c->outside )
           return true; 
        else
           return false; 
     }
     else if( lt == Triangulation::EDGE )
     {
        if( is_outside_VF(triang, Edge(c, u, v)) )
           return true; 
        else
           return false;
     }
     else
     {
        CGAL_assertion( lt == Triangulation::VERTEX );
        return false; 
     }
   }
}

int main()
{
  Triangulation t;

  Vector midpoint(0.5, 0.5);

  Face_handle fh;
  Triangulation::Locate_type lt;
  int i;
  fh = t.locate(Point(0, 0) + midpoint, lt, i);
  CGAL_assertion(lt == Triangulation::EMPTY);

  Vertex_handle vh_midpoint = t.insert(Point(0, 0) + midpoint);
  fh = t.locate(Point(0, 0) + midpoint, lt, i);
  CGAL_assertion(lt == Triangulation::VERTEX && fh->vertex(i) == vh_midpoint);
  t.remove(vh_midpoint);
  CGAL_assertion(t.empty());

  // High degree vertex
  for (int n = 3; n < 8; ++n)
    {
      vh_midpoint = t.insert(Point(0, 0) + midpoint);
      for (int i = 0; i < n; ++i)
        {
          t.insert(Point(0.3 * sin(i * 1.0 / n * 2 * M_PI), 0.3 * cos(i * 1.0 / n * 2 * M_PI)) + midpoint);
        }
      t.remove(vh_midpoint);
      CGAL_assertion(t.is_valid(true));
      while (!t.empty())
        {
          t.remove(t.vertices_begin());
          CGAL_assertion(t.is_valid(true));
        }
    }

  Random random(1284141159);
  std::cout << "Seed: " << random.get_seed () << std::endl;
  Random_points_in_square g(0.495, random);

  CGAL_assertion(t.is_valid());

  std::cout << "Removing first point" << std::endl;
  Vertex_handle vh0 = t.insert(Point(0.5, 0.5));
  CGAL_assertion(t.is_valid());
  t.remove(vh0);
  CGAL_assertion(t.is_valid());
  CGAL_assertion(t.empty());

  {
    Random random(1284141159);
    std::cout << "Seed: " << random.get_seed () << std::endl;

    Random_points_in_square g(0.495, random);
    Vector midpoint(0.5, 0.5);

    Triangulation t;
    CGAL_assertion(t.is_valid());

    std::cout << "Removing first point" << std::endl;
    Vertex_handle vh0 = t.insert(Point(0.5, 0.5));
    CGAL_assertion(t.is_valid());
    t.remove(vh0);
    CGAL_assertion(t.is_valid());
    CGAL_assertion(t.empty());

    std::cout << "Inserting random points and removing them." << std::endl;

    for (int i = 0; i < N_PTS; ++i)
      {
        t.insert(*(++g) + midpoint);
      }
    CGAL_assertion(t.is_valid());

    for (int i = 0; i < N_PTS; ++i)
      {
        // Find a random vertex
        Vertex_handle vh = t.locate(*(++g) + midpoint)->vertex(0);
        vh = t.get_original_vertex(vh);
        t.remove(vh);
        CGAL_assertion(t.is_valid());
      }
  }
  return 0;
}

float
sdf( const Point& q, const Mesh &mesh, const vector<double>& weights,
     KdTree& kd_tree, const Triangulation& triang )
{
  VECTOR3 query (CGAL::to_double(q.x()),
                 CGAL::to_double(q.y()),
                 CGAL::to_double(q.z()));
  kd_tree.queryPosition(query);
  // Initialize the search structure, and search all N points
  int n_vid = kd_tree.getNeighbourPositionIndex(0);
  if(n_vid == -1) throw std::runtime_error("No nearest neighbor. MDS empty?");
  CGAL_assertion( ! mesh.vert_list[n_vid].iso());
  MVertex nv = mesh.vert_list[n_vid];
  double min_sq_d = HUGE;
  int n_fid = -1;
  for(int i = 0; i < nv.num_inc_face; i ++)
  {
     MFace f = mesh.face_list[nv.inc_face(i)];
     Point p[3] = {mesh.vert_list[f.get_corner(0)].point(),
                   mesh.vert_list[f.get_corner(1)].point(),
                   mesh.vert_list[f.get_corner(2)].point()};
     Triangle_3 t (p[0], p[1], p[2]);
     Plane_3 H (p[0], p[1], p[2]);
     Point _q = H.projection(q);
     // check if _q is inside t.
     if( t.has_on(_q) )
     {
        double sq_d = CGAL::to_double((q-_q)*(q-_q));
        if( sq_d < min_sq_d ) { min_sq_d = sq_d; n_fid = nv.inc_face(i); }
     }
     else
     {
        for(int j = 0; j < 3; j ++)
        {
           double _d = CGAL::to_double(CGAL::squared_distance(_q,Segment(p[j], p[(j+1)%3])));
           double sq_d = CGAL::to_double((q-_q)*(q-_q)) + _d;
           if( sq_d < min_sq_d ) { min_sq_d = sq_d; n_fid = nv.inc_face(i); }
        }
     }
  }

  // locate the query point in the triang which is already tagged
  // with in-out flag by the reconstruction.
  bool is_q_outside = false;
  Triangulation::Locate_type lt;
  int u = -1, v = -1;
  Cell_handle c = triang.locate(q, lt, u, v);
  if( lt == Triangulation::OUTSIDE_AFFINE_HULL )
  {
     is_q_outside = true;
     cerr << "Point " << q << " is outside the affine hull." << endl;
  }
  else if( lt == Triangulation::OUTSIDE_CONVEX_HULL )
     is_q_outside = true;
  else
  {
     if( lt == Triangulation::CELL )
     {
        if( c->outside )
           is_q_outside = true;
        else
           is_q_outside = false;
     }
     else if( lt == Triangulation::FACET )
     {
        Cell_handle _c = c->neighbor(u);
        if( c->outside && _c->outside )
           is_q_outside = true;
        else
           is_q_outside = false;
     }
     else if( lt == Triangulation::EDGE )
     {
        if( is_outside_VF(triang, Edge(c, u, v)) )
           is_q_outside = true;
        else
           is_q_outside = false;
     }
     else
     {
        CGAL_assertion( lt == Triangulation::VERTEX );
        is_q_outside = false;
     }
  }

  double w;
  if(weights.size() && mesh.face_list[n_fid].label != -1) w = weights[mesh.face_list[n_fid].label];
  else w = 1.0;
  // double w = mesh.face_list[n_fid].w;

  double gen_sdf = 0;
  if( is_q_outside ) gen_sdf = w*sqrt(min_sq_d);
  else gen_sdf = -w*min_sq_d;

#if 0
  double MAX = 10, MIN = -10;

  if( gen_sdf > MAX ) gen_sdf = MAX;
  if( gen_sdf < MIN ) gen_sdf = MIN;
#endif
//.........这里部分代码省略.........