mirror of https://github.com/CGAL/cgal
Merge branch 'gsoc2013-Visibility_doc-hemmer' of ssh://scm.cgal.org/var/git/cgal-gsoc into gsoc2013-Visibility_doc-hemmer
This commit is contained in:
Michael Hemmer
commit
dedae8798c
|
|
@ -0,0 +1,106 @@
|
|||
// Copyright (c) 2013 Technical University Braunschweig (Germany).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s): Kan Huang <huangkandiy@gmail.com>
|
||||
//
|
||||
|
||||
#ifndef CGAL_PREPROCESSED_VISIBILITY_2_H
|
||||
#define CGAL_PREPROCESSED_VISIBILITY_2_H
|
||||
|
||||
#include <CGAL/Arrangement_2.h>
|
||||
#include <CGAL/Arr_linear_traits_2.h>
|
||||
#include <stack>
|
||||
#include <deque>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template<class Arrangement_2>
|
||||
class Preprocessed_visibility_2 {
|
||||
|
||||
public:
|
||||
typedef typename Arrangement_2::Geometry_traits_2 Geometry_traits_2;
|
||||
// Currently only consider with same type for both
|
||||
typedef Arrangement_2 Input_Arrangement_2;
|
||||
typedef Arrangement_2 Output_Arrangement_2;
|
||||
|
||||
typedef typename Arrangement_2::Halfedge_const_handle Halfedge_const_handle;
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator Ccb_halfedge_const_circulator;
|
||||
typedef typename Arrangement_2::Face_const_handle Face_const_handle;
|
||||
typedef typename Arrangement_2::Kernel Kernel;
|
||||
typedef typename CGAL::Arr_linear_traits_2<Kernel> Linear_traits_2;
|
||||
|
||||
typedef typename Geometry_traits_2::Point_2 Point_2;
|
||||
typedef typename Geometry_traits_2::Ray_2 Ray_2;
|
||||
typedef typename Geometry_traits_2::Segment_2 Segment_2;
|
||||
typedef typename Geometry_traits_2::Line_2 Line_2;
|
||||
typedef typename Geometry_traits_2::Vector_2 Vector_2;
|
||||
typedef typename Geometry_traits_2::FT Number_type;
|
||||
|
||||
typedef typename CGAL::Arrangement_2<Linear_traits_2> Line_Arrangement_2;
|
||||
|
||||
Preprocessed_visibility_2() : p_arr(NULL) {};
|
||||
|
||||
/*! Constructor given an arrangement and the Regularization tag. */
|
||||
Preprocessed_visibility_2(Input_Arrangement_2& arr/*, Regularization_tag r_t*/): p_arr(&arr) {};
|
||||
|
||||
bool is_attached() {
|
||||
return (p_arr != NULL);
|
||||
}
|
||||
|
||||
void attach(Input_Arrangement_2& arr) {
|
||||
p_arr = &arr;
|
||||
}
|
||||
|
||||
void detach() {
|
||||
p_arr = NULL;
|
||||
}
|
||||
|
||||
Input_Arrangement_2 arr() {
|
||||
return *p_arr;
|
||||
}
|
||||
|
||||
void compute_visibility(const Point_2& q,
|
||||
const Face_const_handle face,
|
||||
Output_Arrangement_2& out_arr
|
||||
) {
|
||||
|
||||
}
|
||||
|
||||
void compute_visibility(const Point_2& q,
|
||||
const Halfedge_const_handle he,
|
||||
Output_Arrangement_2& out_arr
|
||||
) {
|
||||
|
||||
}
|
||||
|
||||
private:
|
||||
Input_Arrangement_2* arr;
|
||||
Line_Arrangement_2 line_arr;
|
||||
void preprocess() {
|
||||
|
||||
}
|
||||
|
||||
Line_2 dual_line(const Point_2& p) {
|
||||
return Line_2(p.x(), -1, -p.y());
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif
|
||||
|
|
@ -0,0 +1,866 @@
|
|||
// Copyright (c) 2013 Technical University Braunschweig (Germany).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s): Kan Huang <huangkandiy@gmail.com>
|
||||
//
|
||||
|
||||
#ifndef CGAL_ROTATIONAL_SWEEP_VISIBILITY_2_H
|
||||
#define CGAL_ROTATIONAL_SWEEP_VISIBILITY_2_H
|
||||
|
||||
#include <CGAL/Visibility_2/visibility_utils.h>
|
||||
#include <CGAL/Arrangement_2.h>
|
||||
#include <CGAL/bounding_box.h>
|
||||
#include <boost/unordered_map.hpp>
|
||||
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template<class Arrangement_2_ , class RegularizationTag = CGAL::Tag_true >
|
||||
class Rotational_sweep_visibility_2 {
|
||||
public:
|
||||
typedef Arrangement_2_ Arrangement_2;
|
||||
typedef typename Arrangement_2::Traits_2 Traits_2;
|
||||
typedef typename Arrangement_2::Geometry_traits_2 Geometry_traits_2;
|
||||
typedef typename Arrangement_2::Vertex_const_handle Vertex_const_handle;
|
||||
typedef typename Arrangement_2::Vertex_handle Vertex_handle;
|
||||
typedef typename Arrangement_2::Halfedge_const_handle Halfedge_const_handle;
|
||||
typedef typename Arrangement_2::Halfedge_handle Halfedge_handle;
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator
|
||||
Ccb_halfedge_const_circulator;
|
||||
typedef typename Arrangement_2::Face_const_handle Face_const_handle;
|
||||
typedef typename Arrangement_2::Face_handle Face_handle;
|
||||
|
||||
typedef typename Geometry_traits_2::Kernel K;
|
||||
typedef typename Geometry_traits_2::Point_2 Point_2;
|
||||
typedef typename Geometry_traits_2::Ray_2 Ray_2;
|
||||
typedef typename Geometry_traits_2::Segment_2 Segment_2;
|
||||
typedef typename Geometry_traits_2::Line_2 Line_2;
|
||||
typedef typename Geometry_traits_2::Vector_2 Vector_2;
|
||||
typedef typename Geometry_traits_2::Direction_2 Direction_2;
|
||||
typedef typename Geometry_traits_2::FT Number_type;
|
||||
typedef typename Geometry_traits_2::Object_2 Object_2;
|
||||
|
||||
typedef RegularizationTag Regularization_tag;
|
||||
typedef CGAL::Tag_true Supports_general_polygon_tag;
|
||||
typedef CGAL::Tag_true Supports_simple_polygon_tag;
|
||||
|
||||
private:
|
||||
typedef std::vector<Point_2> Points;
|
||||
typedef Vertex_const_handle VH;
|
||||
typedef std::vector<VH> VHs;
|
||||
typedef Halfedge_const_handle EH;
|
||||
typedef std::vector<EH> EHs;
|
||||
|
||||
class Less_edge: public std::binary_function<EH, EH, bool> {
|
||||
const Geometry_traits_2* geom_traits;
|
||||
public:
|
||||
Less_edge() {}
|
||||
Less_edge(const Geometry_traits_2* traits):geom_traits(traits) {}
|
||||
bool operator() (const EH e1, const EH e2) const {
|
||||
if (e1 == e2)
|
||||
return false;
|
||||
else {
|
||||
return &(*e1)<&(*e2);
|
||||
// if (e1->source() == e2->source())
|
||||
// return Visibility_2::compare_xy_2(geom_traits, e1->target()->point(), e2->target()->point()) == SMALLER;
|
||||
// else
|
||||
// return Visibility_2::compare_xy_2(geom_traits, e1->source()->point(), e2->source()->point()) == SMALLER;
|
||||
}
|
||||
}
|
||||
};
|
||||
|
||||
class Less_vertex: public std::binary_function<VH, VH, bool> {
|
||||
const Geometry_traits_2* geom_traits;
|
||||
public:
|
||||
Less_vertex() {}
|
||||
Less_vertex(const Geometry_traits_2* traits):geom_traits(traits) {}
|
||||
bool operator() (const VH v1, const VH v2) const {
|
||||
if (v1 == v2)
|
||||
return false;
|
||||
else
|
||||
// I know this is dirty but it speeds up by 25%. Michael
|
||||
return &(*v1)<&(*v2);
|
||||
// return Visibility_2::compare_xy_2(geom_traits, v1->point(), v2->point()) == SMALLER;
|
||||
}
|
||||
};
|
||||
|
||||
class Closer_edge: public std::binary_function<EH, EH, bool> {
|
||||
const Geometry_traits_2* geom_traits;
|
||||
Point_2 q;
|
||||
public:
|
||||
Closer_edge() {}
|
||||
Closer_edge(const Geometry_traits_2* traits, const Point_2& q):geom_traits(traits), q(q) {}
|
||||
|
||||
int vtype(const Point_2& c, const Point_2& p) const {
|
||||
switch(Visibility_2::orientation_2(geom_traits, q, c, p)) {
|
||||
case COLLINEAR:
|
||||
if (Visibility_2::less_distance_to_point_2(geom_traits, q, c, p))
|
||||
return 0;
|
||||
else
|
||||
return 3;
|
||||
case RIGHT_TURN:
|
||||
return 1;
|
||||
case LEFT_TURN:
|
||||
return 2;
|
||||
}
|
||||
}
|
||||
bool operator() (const EH& e1, const EH& e2) const {
|
||||
if (e1 == e2)
|
||||
return false;
|
||||
const Point_2& s1=e1->source()->point(),
|
||||
t1=e1->target()->point(),
|
||||
s2=e2->source()->point(),
|
||||
t2=e2->target()->point();
|
||||
if (e1->source() == e2->source()) {
|
||||
|
||||
int vt1 = vtype(s1, t1),
|
||||
vt2 = vtype(s1, t2);
|
||||
if (vt1 != vt2)
|
||||
return vt1 > vt2;
|
||||
else
|
||||
return (Visibility_2::orientation_2(geom_traits, s1, t2, t1)==
|
||||
Visibility_2::orientation_2(geom_traits, s1, t2, q));
|
||||
}
|
||||
|
||||
if (e1->target() == e2->source()) {
|
||||
// const Point_2& p1 = s1,
|
||||
// p2 = t2,
|
||||
// c = s2;
|
||||
int vt1 = vtype(t1, s1),
|
||||
vt2 = vtype(t1, t2);
|
||||
if (vt1 != vt2)
|
||||
return vt1 > vt2;
|
||||
else
|
||||
return (Visibility_2::orientation_2(geom_traits, s2, t2, s1)==
|
||||
Visibility_2::orientation_2(geom_traits, s2, t2, q));
|
||||
}
|
||||
|
||||
if (e1->source() == e2->target()) {
|
||||
// const Point_2& p1 = t1,
|
||||
// p2 = s2,
|
||||
// c = s1;
|
||||
int vt1 = vtype(s1, t1),
|
||||
vt2 = vtype(s1, s2);
|
||||
if (vt1 != vt2)
|
||||
return vt1 > vt2;
|
||||
else return (Visibility_2::orientation_2(geom_traits, s1, s2, t1)==
|
||||
Visibility_2::orientation_2(geom_traits, s1, s2, q));
|
||||
}
|
||||
|
||||
if (e1->target() == e2->target()) {
|
||||
// const Point_2& p1 = s1,
|
||||
// p2 = s2,
|
||||
// c = t1;
|
||||
int vt1 = vtype(t1, s1),
|
||||
vt2 = vtype(t1, s2);
|
||||
if (vt1 != vt2)
|
||||
return vt1 > vt2;
|
||||
else return (Visibility_2::orientation_2(geom_traits, t1, s2, s1)==
|
||||
Visibility_2::orientation_2(geom_traits, t1, s2, q));
|
||||
}
|
||||
|
||||
Orientation e1q = Visibility_2::orientation_2(geom_traits, s1, t1, q);
|
||||
switch (e1q)
|
||||
{
|
||||
case COLLINEAR:
|
||||
if (Visibility_2::collinear(geom_traits, q, s2, t2)) {
|
||||
//q is collinear with e1 and e2.
|
||||
return (Visibility_2::less_distance_to_point_2(geom_traits, q, s1, s2)
|
||||
|| Visibility_2::less_distance_to_point_2(geom_traits, q, t1, t2));
|
||||
}
|
||||
else {
|
||||
//q is collinear with e1 not with e2.
|
||||
if (Visibility_2::collinear(geom_traits, s2, t2, s1))
|
||||
return (Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, t1));
|
||||
else
|
||||
return (Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, s1));
|
||||
}
|
||||
case RIGHT_TURN:
|
||||
switch (Visibility_2::orientation_2(geom_traits, s1, t1, s2)) {
|
||||
case COLLINEAR:
|
||||
return Visibility_2::orientation_2(geom_traits, s1, t1, t2)!=e1q;
|
||||
case RIGHT_TURN:
|
||||
if (Visibility_2::orientation_2(geom_traits, s1, t1, t2) == LEFT_TURN)
|
||||
return Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, s1);
|
||||
else
|
||||
return false;
|
||||
case LEFT_TURN:
|
||||
if (Visibility_2::orientation_2(geom_traits, s1, t1, t2) == RIGHT_TURN)
|
||||
return Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, s1);
|
||||
else
|
||||
return true;
|
||||
}
|
||||
case LEFT_TURN:
|
||||
switch (Visibility_2::orientation_2(geom_traits, s1, t1, s2)) {
|
||||
case COLLINEAR:
|
||||
return Visibility_2::orientation_2(geom_traits, s1, t1, t2)!=e1q;
|
||||
case LEFT_TURN:
|
||||
if (Visibility_2::orientation_2(geom_traits, s1, t1, t2) == RIGHT_TURN)
|
||||
return Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, s1);
|
||||
else
|
||||
return false;
|
||||
case RIGHT_TURN:
|
||||
if (Visibility_2::orientation_2(geom_traits, s1, t1, t2) == LEFT_TURN)
|
||||
return Visibility_2::orientation_2(geom_traits, s2, t2, q)
|
||||
== Visibility_2::orientation_2(geom_traits, s2, t2, s1);
|
||||
else
|
||||
return true;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
// Using hash_map or edx causes a seg fault, did not have the time to see why. Michael
|
||||
// class Hash_edge: public std::unary_function<VH,typename boost::hash<const typename Arrangement_2::X_monotone_curve_2*>::result_type> {
|
||||
// public:
|
||||
// typename boost::hash<const typename Arrangement_2::X_monotone_curve_2*>::result_type
|
||||
// operator() (const EH e1) const {
|
||||
// return boost::hash<const typename Arrangement_2::X_monotone_curve_2*>()(&(e1->curve()));
|
||||
// }
|
||||
// };
|
||||
|
||||
|
||||
const Geometry_traits_2 *geom_traits;
|
||||
const Arrangement_2 *p_arr;
|
||||
Point_2 q; //query point
|
||||
Points polygon; //visibility polygon
|
||||
std::map<VH, EHs, Less_vertex> incident_edges; //the edges that are
|
||||
std::map<EH, int, Less_edge> edx; //index of active edges in the heap
|
||||
// boost::unordered_map<EH,int,Hash_edge> edx; //index of active edges in the heap
|
||||
std::set<EH, Closer_edge> active_edges; //a set of edges that intersect the current vision ray.
|
||||
VHs vs; //angular sorted vertices
|
||||
EHs bad_edges; //edges that pass the query point
|
||||
VH cone_end1; //an end of visibility cone
|
||||
VH cone_end2; //another end of visibility cone
|
||||
int cone_end1_idx; //index of cone_end1->point() in visibility polygon
|
||||
int cone_end2_idx; //index of cone_end2->point() in visibility polygon
|
||||
|
||||
bool is_vertex_query;
|
||||
bool is_edge_query;
|
||||
bool is_face_query;
|
||||
bool is_big_cone; //whether the angle of visibility_cone is greater than pi.
|
||||
|
||||
public:
|
||||
Rotational_sweep_visibility_2(): p_arr(NULL), geom_traits(NULL) {}
|
||||
Rotational_sweep_visibility_2(const Arrangement_2& arr): p_arr(&arr) {
|
||||
geom_traits = p_arr->geometry_traits();
|
||||
}
|
||||
|
||||
const std::string name(){return std::string("R_visibility_2");}
|
||||
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(const Point_2& q, const Halfedge_const_handle e, VARR& arr_out) {
|
||||
arr_out.clear();
|
||||
bad_edges.clear();
|
||||
this->q = q;
|
||||
|
||||
if (Visibility_2::compare_xy_2(geom_traits, q, e->target()->point())==EQUAL) {
|
||||
is_vertex_query = true;
|
||||
is_edge_query = false;
|
||||
is_face_query = false;
|
||||
cone_end1 = e->source();
|
||||
cone_end2 = e->next()->target();
|
||||
is_big_cone = CGAL::right_turn(cone_end1->point(), q, cone_end2->point());
|
||||
|
||||
typename Arrangement_2::Halfedge_around_vertex_const_circulator first, curr;
|
||||
first = curr = e->target()->incident_halfedges();
|
||||
do {
|
||||
if (curr->face() == e->face())
|
||||
bad_edges.push_back(curr);
|
||||
else if (curr->twin()->face() == e->face())
|
||||
bad_edges.push_back(curr->twin());
|
||||
} while (++curr != first);
|
||||
}
|
||||
else {
|
||||
is_vertex_query = false;
|
||||
is_edge_query = true;
|
||||
is_face_query = false;
|
||||
cone_end1 = e->source();
|
||||
cone_end2 = e->target();
|
||||
bad_edges.push_back(e);
|
||||
is_big_cone = false;
|
||||
}
|
||||
visibility_region_impl(e->face(), q);
|
||||
|
||||
//decide which inside of the visibility butterfly is needed.
|
||||
int small_idx, big_idx;
|
||||
if ( cone_end1_idx < cone_end2_idx ) {
|
||||
small_idx = cone_end1_idx;
|
||||
big_idx = cone_end2_idx;
|
||||
}
|
||||
else {
|
||||
small_idx = cone_end2_idx;
|
||||
big_idx = cone_end1_idx;
|
||||
}
|
||||
int next_idx = small_idx + 1;
|
||||
bool is_between;
|
||||
//indicate whether the shape between small_idx and big_idx is the visibility region required.
|
||||
if (CGAL::right_turn(cone_end1->point(), q, cone_end2->point())) {
|
||||
is_between = false;
|
||||
while (next_idx != big_idx) {
|
||||
if (CGAL::left_turn(cone_end1->point(), q, polygon[next_idx]) || CGAL::left_turn(q, cone_end2->point(), polygon[next_idx])) {
|
||||
is_between = true;
|
||||
break;
|
||||
}
|
||||
next_idx++;
|
||||
}
|
||||
}
|
||||
else {
|
||||
is_between = true;
|
||||
while (next_idx != big_idx) {
|
||||
if (CGAL::right_turn(cone_end1->point(), q, polygon[next_idx]) || CGAL::right_turn(q, cone_end2->point(), polygon[next_idx])) {
|
||||
is_between = false;
|
||||
break;
|
||||
}
|
||||
next_idx++;
|
||||
}
|
||||
}
|
||||
|
||||
typename Points::iterator first = polygon.begin() + small_idx;
|
||||
typename Points::iterator last = polygon.begin() + big_idx;
|
||||
if (is_between) {
|
||||
Points polygon_out(first, last+1);
|
||||
if (is_vertex_query)
|
||||
polygon_out.push_back(q);
|
||||
Visibility_2::report_while_handling_needles<Rotational_sweep_visibility_2>
|
||||
(geom_traits, q, polygon_out, arr_out);
|
||||
}
|
||||
else {
|
||||
Points polygon_out(polygon.begin(), first+1);
|
||||
if (is_vertex_query) polygon_out.push_back(q);
|
||||
for (int i = big_idx; i != polygon.size(); i++) {
|
||||
polygon_out.push_back(polygon[i]);
|
||||
}
|
||||
Visibility_2::report_while_handling_needles<Rotational_sweep_visibility_2>
|
||||
(geom_traits, q, polygon_out, arr_out);
|
||||
}
|
||||
|
||||
conditional_regularize(arr_out, Regularization_tag());
|
||||
|
||||
if (arr_out.faces_begin()->is_unbounded())
|
||||
return ++arr_out.faces_begin();
|
||||
else
|
||||
return arr_out.faces_begin();
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(const Point_2& q, const Face_const_handle f, VARR& arr_out) {
|
||||
arr_out.clear();
|
||||
this->q = q;
|
||||
is_vertex_query = false;
|
||||
is_edge_query = false;
|
||||
is_face_query = true;
|
||||
|
||||
visibility_region_impl(f, q);
|
||||
Visibility_2::report_while_handling_needles<Rotational_sweep_visibility_2>(geom_traits, q, polygon, arr_out);
|
||||
conditional_regularize(arr_out, Regularization_tag());
|
||||
if (arr_out.faces_begin()->is_unbounded())
|
||||
return ++arr_out.faces_begin();
|
||||
else
|
||||
return arr_out.faces_begin();
|
||||
}
|
||||
|
||||
bool is_attached() {
|
||||
return (p_arr != NULL);
|
||||
}
|
||||
|
||||
void attach(const Arrangement_2& arr) {
|
||||
p_arr = &arr;
|
||||
geom_traits = p_arr->geometry_traits();
|
||||
}
|
||||
|
||||
void detach() {
|
||||
p_arr = NULL;
|
||||
geom_traits = NULL;
|
||||
}
|
||||
|
||||
const Arrangement_2& arr() {
|
||||
return *p_arr;
|
||||
}
|
||||
|
||||
private:
|
||||
//get the neighbor of v along edge e
|
||||
VH get_neighbor(const EH e, const VH v) {
|
||||
if (e->source() == v)
|
||||
return e->target();
|
||||
else
|
||||
return e->source();
|
||||
}
|
||||
//check whether ray(q->dp) intersects segment(p1, p2)
|
||||
bool do_intersect_ray(const Point_2& q,
|
||||
const Point_2& dp,
|
||||
const Point_2& p1,
|
||||
const Point_2& p2) {
|
||||
return (CGAL::orientation(q, dp, p1) != CGAL::orientation(q, dp, p2) && CGAL::orientation(q, p1, dp) == CGAL::orientation(q, p1, p2));
|
||||
}
|
||||
|
||||
//arrange vertices that on a same vision ray in a 'funnel' order
|
||||
void funnel(int i, int j) {
|
||||
VHs right, left;
|
||||
//whether the edges incident to a vertex block the left side and right side of current vision ray.
|
||||
bool block_left(false), block_right(false);
|
||||
VH former = vs[i], nb;
|
||||
for (int l=i; l<j; l++) {
|
||||
bool left_v(false), right_v(false), has_predecessor(false);
|
||||
EHs& edges = incident_edges[vs[l]];
|
||||
for (int k=0; k<edges.size(); k++) {
|
||||
nb = get_neighbor(edges[k], vs[l]);
|
||||
if ( nb == former ) {
|
||||
has_predecessor = true;
|
||||
break;
|
||||
}
|
||||
if (CGAL::left_turn(q, vs[l]->point(), nb->point()))
|
||||
left_v = true;
|
||||
else
|
||||
right_v = CGAL::right_turn(q, vs[l]->point(), nb->point());
|
||||
}
|
||||
if (has_predecessor) {
|
||||
//if the current vertex connects to the vertex before by an edge,
|
||||
//the vertex before can help it to block.
|
||||
block_left = block_left || left_v;
|
||||
block_right = block_right || right_v;
|
||||
}
|
||||
else {
|
||||
block_left = left_v;
|
||||
block_right = right_v;
|
||||
}
|
||||
if (block_left && block_right) {
|
||||
//when both sides are blocked, there is no need to change the vertex after.
|
||||
right.push_back(vs[l]);
|
||||
break;
|
||||
}
|
||||
else {
|
||||
if (block_left)
|
||||
left.push_back(vs[l]);
|
||||
else
|
||||
right.push_back(vs[l]);
|
||||
}
|
||||
former = vs[l];
|
||||
}
|
||||
for (int l=0; l!=right.size(); l++)
|
||||
vs[i+l] = right[l];
|
||||
for (int l=0; l!=left.size(); l++)
|
||||
vs[i+l+right.size()] = left[left.size()-1-l];
|
||||
}
|
||||
|
||||
|
||||
|
||||
void visibility_region_impl(const Face_const_handle f, const Point_2& q) {
|
||||
vs.clear();
|
||||
polygon.clear();
|
||||
active_edges = std::set<EH, Closer_edge>(Closer_edge(geom_traits, q));
|
||||
incident_edges = std::map<VH, EHs, Less_vertex>(Less_vertex(geom_traits));
|
||||
edx = std::map<EH, int, Less_edge>(Less_edge(geom_traits));
|
||||
|
||||
EHs relevant_edges; //all edges that can affect the visibility of query point.
|
||||
Arrangement_2 bbox;
|
||||
if (is_face_query)
|
||||
input_face(f);
|
||||
else
|
||||
input_face(f, relevant_edges, bbox);
|
||||
//the following code is the initiation of vision ray. the direction of the initial ray is between the direction
|
||||
//from q to last vertex in vs and positive x-axis. By choosing this direction, we make
|
||||
//sure that all plane is swept and there is not needle at the beginning of sweeping.
|
||||
Vector_2 dir;
|
||||
if (Direction_2(-1, 0) < Direction_2(Vector_2(q, vs.back()->point())))
|
||||
dir = Vector_2(1, 0) + Vector_2(q, vs.back()->point());
|
||||
else
|
||||
dir = Vector_2(0, -1);
|
||||
Point_2 dp = q + dir;
|
||||
|
||||
//initiation of active_edges. for face queries, all edges on the boundary can affect visibility.
|
||||
//for non-face queries, only relevant_edges has to be considered.
|
||||
if (is_face_query) {
|
||||
Ccb_halfedge_const_circulator curr = f->outer_ccb();
|
||||
Ccb_halfedge_const_circulator circ = curr;
|
||||
do {
|
||||
if (do_intersect_ray(q, dp, curr->target()->point(), curr->source()->point())) {
|
||||
active_edges.insert(curr);
|
||||
}
|
||||
} while (++curr != circ);
|
||||
|
||||
typename Arrangement_2::Hole_const_iterator hi;
|
||||
for (hi = f->holes_begin(); hi != f->holes_end(); ++hi) {
|
||||
Ccb_halfedge_const_circulator curr = circ = *hi;
|
||||
do {
|
||||
if (do_intersect_ray(q, dp, curr->target()->point(), curr->source()->point()))
|
||||
active_edges.insert(curr);
|
||||
} while (++curr != circ);
|
||||
}
|
||||
}
|
||||
else {
|
||||
for (int i=0; i!=relevant_edges.size(); i++)
|
||||
if (do_intersect_ray(q, dp, relevant_edges[i]->source()->point(), relevant_edges[i]->target()->point()))
|
||||
active_edges.insert(relevant_edges[i]);
|
||||
}
|
||||
|
||||
//angular sweep begins
|
||||
// std::cout<<active_edges.size()<<std::endl;
|
||||
for (int i=0; i!=vs.size(); i++) {
|
||||
VH vh = vs[i];
|
||||
EH closest_e = *active_edges.begin();
|
||||
EHs& edges = incident_edges[vh];
|
||||
EHs insert_ehs, remove_ehs;
|
||||
for (int j=0; j!=edges.size(); j++) {
|
||||
EH& e = edges[j];
|
||||
if (active_edges.find(e) == active_edges.end())
|
||||
insert_ehs.push_back(e);
|
||||
else
|
||||
remove_ehs.push_back(e);
|
||||
}
|
||||
int insert_cnt = insert_ehs.size();
|
||||
int remove_cnt = remove_ehs.size();
|
||||
if (insert_cnt == 1 && remove_cnt == 1) {
|
||||
const EH& ctemp_eh = *active_edges.find(remove_ehs.front());
|
||||
EH& temp_eh = const_cast<EH&>(ctemp_eh);
|
||||
temp_eh = insert_ehs.front();
|
||||
}
|
||||
else {
|
||||
for (int j=0; j!=remove_cnt; j++)
|
||||
active_edges.erase(remove_ehs[j]);
|
||||
for (int j=0; j!=insert_cnt; j++)
|
||||
active_edges.insert(insert_ehs[j]);
|
||||
}
|
||||
|
||||
if (closest_e != *active_edges.begin()) {
|
||||
//when the closest edge changed
|
||||
if (is_face_query) {
|
||||
if (remove_cnt > 0 && insert_cnt > 0) {
|
||||
//some edges are added and some are deleted, which means the vertex swept is part of visibility polygon.
|
||||
update_visibility(vh->point());
|
||||
}
|
||||
if (remove_cnt == 0 && insert_cnt > 0) {
|
||||
//only add some edges, means the view ray is blocked by new edges.
|
||||
//therefore first add the intersection of view ray and former closet edge, then add the vertice swept.
|
||||
update_visibility(ray_seg_intersection(q,
|
||||
vh->point(),
|
||||
closest_e->target()->point(),
|
||||
closest_e->source()->point()));
|
||||
update_visibility(vh->point());
|
||||
}
|
||||
if (remove_cnt > 0 && insert_cnt == 0) {
|
||||
//only delete some edges, means some block is moved and the view ray can reach the segments after the block.
|
||||
update_visibility(vh->point());
|
||||
update_visibility(ray_seg_intersection(q,
|
||||
vh->point(),
|
||||
(*active_edges.begin())->target()->point(),
|
||||
(*active_edges.begin())->source()->point()));
|
||||
}
|
||||
}
|
||||
else {
|
||||
//extra work here for edge/vertex query is the index of cone_end1 and cone_end2 will be recorded.
|
||||
if (remove_cnt > 0 && insert_cnt > 0) {
|
||||
//some edges are added and some are deleted, which means the vertice swept is part of visibility polygon.
|
||||
if (update_visibility(vh->point())) {
|
||||
if (vh == cone_end1)
|
||||
cone_end1_idx = polygon.size()-1;
|
||||
else if (vh == cone_end2)
|
||||
cone_end2_idx = polygon.size()-1;
|
||||
}
|
||||
}
|
||||
if (remove_cnt == 0 && insert_cnt > 0) {
|
||||
//only add some edges, means the view ray is blocked by new edges.
|
||||
//therefore first add the intersection of view ray and former closet edge, then add the vertice swept.
|
||||
update_visibility(ray_seg_intersection(q,
|
||||
vh->point(),
|
||||
closest_e->target()->point(),
|
||||
closest_e->source()->point()));
|
||||
if (update_visibility(vh->point())) {
|
||||
if (vh == cone_end1)
|
||||
cone_end1_idx = polygon.size()-1;
|
||||
else if (vh == cone_end2)
|
||||
cone_end2_idx = polygon.size()-1;
|
||||
}
|
||||
}
|
||||
if (remove_cnt > 0 && insert_cnt == 0) {
|
||||
//only delete some edges, means some block is removed and the vision ray can reach the segments after the block.
|
||||
if (update_visibility(vh->point())) {
|
||||
if (vh == cone_end1)
|
||||
cone_end1_idx = polygon.size()-1;
|
||||
else if (vh == cone_end2)
|
||||
cone_end2_idx = polygon.size()-1;
|
||||
}
|
||||
update_visibility(ray_seg_intersection(q,
|
||||
vh->point(),
|
||||
(*active_edges.begin())->target()->point(),
|
||||
(*active_edges.begin())->source()->point()));
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void print_edge(const EH e) {
|
||||
std::cout<<e->source()->point()<<"->"<<e->target()->point()<<std::endl;
|
||||
}
|
||||
|
||||
//compute the intersection of ray(q->dp) and segment(s, t)
|
||||
//if they are collinear then return the endpoint which is closer to q.
|
||||
|
||||
Point_2 ray_seg_intersection(
|
||||
const Point_2& q, const Point_2& dp, // the ray
|
||||
const Point_2& s, const Point_2& t) // the segment
|
||||
{
|
||||
if (CGAL::collinear(q, dp, s)) {
|
||||
if (CGAL::collinear(q, dp, t)) {
|
||||
if (CGAL::compare_distance_to_point(q, s, t)==CGAL::SMALLER)
|
||||
return s;
|
||||
else
|
||||
return t;
|
||||
}
|
||||
else
|
||||
return s;
|
||||
}
|
||||
Ray_2 ray(q,dp);
|
||||
Segment_2 seg(s,t);
|
||||
CGAL::Object result = CGAL::intersection(ray, seg);
|
||||
return *(CGAL::object_cast<Point_2>(&result));
|
||||
}
|
||||
|
||||
//check if p has been discovered before, if not update the visibility polygon
|
||||
bool update_visibility(const Point_2& p){
|
||||
if (polygon.empty()) {
|
||||
polygon.push_back(p);
|
||||
return true;
|
||||
}
|
||||
else if (Visibility_2::compare_xy_2(geom_traits, polygon.back(), p) != EQUAL) {
|
||||
polygon.push_back(p);
|
||||
return true;
|
||||
}
|
||||
return false;
|
||||
}
|
||||
|
||||
//functor to decide which vertex is swept earlier by the rotational sweeping ray
|
||||
class Is_swept_earlier:public std::binary_function<VH, VH, bool> {
|
||||
const Point_2& q;
|
||||
const Geometry_traits_2* geom_traits;
|
||||
public:
|
||||
Is_swept_earlier(const Point_2& q, const Geometry_traits_2* traits):q(q), geom_traits(traits) {}
|
||||
bool operator() (const VH v1, const VH v2) const {
|
||||
const Point_2& p1 = v1->point();
|
||||
const Point_2& p2 = v2->point();
|
||||
int qua1 = quadrant(q, p1);
|
||||
int qua2 = quadrant(q, p2);
|
||||
if (qua1 < qua2)
|
||||
return true;
|
||||
if (qua1 > qua2)
|
||||
return false;
|
||||
if (collinear(q, p1, p2))
|
||||
return (CGAL::compare_distance_to_point(q, p1, p2) == CGAL::SMALLER);
|
||||
else
|
||||
return CGAL::right_turn(p1, q, p2);
|
||||
}
|
||||
|
||||
//return the quadrant of p with respect to o.
|
||||
int quadrant(const Point_2& o, const Point_2& p) const {
|
||||
typename Geometry_traits_2::Compare_x_2 compare_x = geom_traits->compare_x_2_object();
|
||||
typename Geometry_traits_2::Compare_y_2 compare_y = geom_traits->compare_y_2_object();
|
||||
|
||||
Comparison_result dx = compare_x(p, o);
|
||||
Comparison_result dy = compare_y(p, o);
|
||||
if (dx==LARGER && dy!=SMALLER)
|
||||
return 1;
|
||||
if (dx!=LARGER && dy==LARGER)
|
||||
return 2;
|
||||
if (dx==SMALLER && dy!=LARGER)
|
||||
return 3;
|
||||
if (dx!=SMALLER && dy==SMALLER)
|
||||
return 4;
|
||||
return 0;
|
||||
}
|
||||
};
|
||||
|
||||
//when the query point is in face, every edge is good.
|
||||
void input_neighbor_f( const Halfedge_const_handle e) {
|
||||
VH v = e->target();
|
||||
if (!incident_edges.count(v))
|
||||
vs.push_back(v);
|
||||
incident_edges[v].push_back(e);
|
||||
incident_edges[v].push_back(e->next());
|
||||
}
|
||||
|
||||
//check if p is in the visibility cone
|
||||
bool is_in_cone(const Point_2& p) const{
|
||||
if (is_big_cone)
|
||||
return (!CGAL::right_turn(cone_end1->point(), q, p)) || (!CGAL::left_turn(cone_end2->point(), q, p));
|
||||
else
|
||||
return (!CGAL::right_turn(cone_end1->point(), q, p)) && (!CGAL::left_turn(cone_end2->point(), q, p));
|
||||
}
|
||||
|
||||
//for vertex and edge query: the visibility is limited in a cone.
|
||||
void input_edge(const Halfedge_const_handle e,
|
||||
EHs& good_edges) {
|
||||
for (int i=0; i<bad_edges.size(); i++)
|
||||
if (e == bad_edges[i])
|
||||
return;
|
||||
VH v1 = e->target();
|
||||
VH v2 = e->source();
|
||||
//an edge will affect visibility only if it has an endpoint in the visibility cone or it crosses the boundary of the cone.
|
||||
if (is_in_cone(v1->point()) || is_in_cone(v2->point()) || do_intersect_ray(q, cone_end1->point(), v1->point(), v2->point())) {
|
||||
good_edges.push_back(e);
|
||||
if (!incident_edges.count(v1))
|
||||
vs.push_back(v1);
|
||||
incident_edges[v1].push_back(e);
|
||||
if (!incident_edges.count(v2))
|
||||
vs.push_back(v2);
|
||||
incident_edges[v2].push_back(e);
|
||||
}
|
||||
}
|
||||
|
||||
//for face query: traverse the face to get all edges
|
||||
//and sort vertices in counter-clockwise order.
|
||||
void input_face (Face_const_handle fh)
|
||||
{
|
||||
Ccb_halfedge_const_circulator curr = fh->outer_ccb();
|
||||
Ccb_halfedge_const_circulator circ = curr;
|
||||
do {
|
||||
assert(curr->face() == fh);
|
||||
input_neighbor_f(curr);
|
||||
} while (++curr != circ);
|
||||
|
||||
typename Arrangement_2::Hole_const_iterator hi;
|
||||
for (hi = fh->holes_begin(); hi != fh->holes_end(); ++hi) {
|
||||
Ccb_halfedge_const_circulator curr = *hi, circ = *hi;
|
||||
do {
|
||||
assert(curr->face() == fh);
|
||||
input_neighbor_f(curr);
|
||||
} while (++curr != circ);
|
||||
}
|
||||
|
||||
std::sort(vs.begin(), vs.end(), Is_swept_earlier(q, geom_traits));
|
||||
|
||||
for (int i=0; i!=vs.size(); i++) {
|
||||
int j = i+1;
|
||||
while (j != vs.size()) {
|
||||
if (!CGAL::collinear(q, vs[i]->point(), vs[j]->point()))
|
||||
break;
|
||||
j++;
|
||||
}
|
||||
if (j-i>1)
|
||||
funnel(i, j);
|
||||
i = j-1;
|
||||
}
|
||||
}
|
||||
//for vertex or edge query: traverse the face to get all edges
|
||||
//and sort vertices in counter-clockwise order.
|
||||
void input_face (Face_const_handle fh,
|
||||
EHs& good_edges,
|
||||
Arrangement_2& bbox)
|
||||
{
|
||||
Ccb_halfedge_const_circulator curr = fh->outer_ccb();
|
||||
Ccb_halfedge_const_circulator circ = curr;
|
||||
do {
|
||||
assert(curr->face() == fh);
|
||||
input_edge(curr, good_edges);
|
||||
} while (++curr != circ);
|
||||
|
||||
typename Arrangement_2::Hole_const_iterator hi;
|
||||
for (hi = fh->holes_begin(); hi != fh->holes_end(); ++hi) {
|
||||
Ccb_halfedge_const_circulator curr = circ = *hi;
|
||||
do {
|
||||
assert(curr->face() == fh);
|
||||
input_edge(curr, good_edges);
|
||||
} while (++curr != circ);
|
||||
}
|
||||
|
||||
//create a box that cover all vertices such that during the sweeping, the vision ray will always intersect at least an edge.
|
||||
//this box doesn't intersect any relevant_edge.
|
||||
Points points;
|
||||
for (int i=0; i<vs.size(); i++) {
|
||||
points.push_back(vs[i]->point());
|
||||
}
|
||||
points.push_back(q);
|
||||
//first get the bounding box of all relevant points.
|
||||
typename Geometry_traits_2::Iso_rectangle_2 bb = bounding_box(points.begin(), points.end());
|
||||
|
||||
Number_type xmin, xmax, ymin, ymax;
|
||||
typename Geometry_traits_2::Compute_x_2 compute_x = geom_traits->compute_x_2_object();
|
||||
typename Geometry_traits_2::Compute_y_2 compute_y = geom_traits->compute_y_2_object();
|
||||
|
||||
//make the box a little bigger than bb so that it won't intersect any relevant_edge.
|
||||
xmin = compute_x(bb.min())-1;
|
||||
ymin = compute_y(bb.min())-1;
|
||||
xmax = compute_x(bb.max())+1;
|
||||
ymax = compute_y(bb.max())+1;
|
||||
Point_2 box[4] = {Point_2(xmin, ymin), Point_2(xmax, ymin),
|
||||
Point_2(xmax, ymax), Point_2(xmin, ymax)};
|
||||
Halfedge_handle e1 = bbox.insert_in_face_interior(Segment_2(box[0], box[1]), bbox.unbounded_face());
|
||||
Halfedge_handle e2 = bbox.insert_from_left_vertex(Segment_2(box[1], box[2]), e1->target());
|
||||
Halfedge_handle e3 = bbox.insert_from_right_vertex(Segment_2(box[2], box[3]), e2->target());
|
||||
bbox.insert_at_vertices(Segment_2(box[0], box[3]), e1->source(), e3->target());
|
||||
|
||||
circ = curr = e1->face()->outer_ccb();
|
||||
do {
|
||||
VH v = curr->target();
|
||||
vs.push_back(v);
|
||||
incident_edges[v].push_back(curr);
|
||||
incident_edges[v].push_back(curr->next());
|
||||
good_edges.push_back(curr);
|
||||
} while(++curr != circ);
|
||||
|
||||
std::sort(vs.begin(), vs.end(), Is_swept_earlier(q, geom_traits));
|
||||
|
||||
for (int i=0; i!=vs.size(); i++) {
|
||||
int j = i+1;
|
||||
while (j != vs.size()) {
|
||||
if (!CGAL::collinear(q, vs[i]->point(), vs[j]->point()))
|
||||
break;
|
||||
j++;
|
||||
}
|
||||
if (j-i>1)
|
||||
funnel(i, j);
|
||||
i = j-1;
|
||||
}
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
void conditional_regularize(VARR& arr_out, CGAL::Tag_true) {
|
||||
regularize_output(arr_out);
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
void conditional_regularize(VARR& arr_out, CGAL::Tag_false) {
|
||||
//do nothing
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
void regularize_output(VARR& arr_out) {
|
||||
typename VARR::Edge_iterator e_itr;
|
||||
for (e_itr = arr_out.edges_begin();
|
||||
e_itr != arr_out.edges_end();
|
||||
e_itr++) {
|
||||
if (e_itr->face() == e_itr->twin()->face())
|
||||
arr_out.remove_edge(e_itr);
|
||||
}
|
||||
}
|
||||
};
|
||||
} // end namespace CGAL
|
||||
|
||||
|
||||
|
||||
#endif
|
||||
|
||||
|
||||
|
|
@ -0,0 +1,715 @@
|
|||
// Copyright (c) 2013 Technical University Braunschweig (Germany).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$`
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s): Francisc Bungiu <fbungiu@gmail.com>
|
||||
// Michael Hemmer <michael.hemmer@cgal.org>
|
||||
// Ning Xu <longyin0904@gmail.com>
|
||||
|
||||
#ifndef CGAL_SIMPLE_POLYGON_VISIBILITY_2_H
|
||||
#define CGAL_SIMPLE_POLYGON_VISIBILITY_2_H
|
||||
|
||||
#include <CGAL/Arrangement_2.h>
|
||||
#include <CGAL/tags.h>
|
||||
#include <CGAL/enum.h>
|
||||
#include <CGAL/Visibility_2/visibility_utils.h>
|
||||
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
|
||||
#include <stack>
|
||||
#include <map>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template<class Arrangement_2_, class RegularizationTag = CGAL::Tag_true>
|
||||
class Simple_polygon_visibility_2 {
|
||||
|
||||
public:
|
||||
// Currently only consider with same type for both
|
||||
typedef Arrangement_2_ Arrangement_2;
|
||||
typedef typename Arrangement_2::Traits_2 Traits_2;
|
||||
typedef typename Arrangement_2::Geometry_traits_2 Geometry_traits_2;
|
||||
typedef typename Geometry_traits_2::Kernel K;
|
||||
|
||||
typedef typename Arrangement_2::Halfedge_const_handle
|
||||
Halfedge_const_handle;
|
||||
typedef typename Arrangement_2::Halfedge_handle Halfedge_handle;
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator
|
||||
Ccb_halfedge_const_circulator;
|
||||
typedef typename Arrangement_2::Face_const_handle Face_const_handle;
|
||||
typedef typename Arrangement_2::Face_handle Face_handle;
|
||||
typedef typename Arrangement_2::Halfedge_around_vertex_const_circulator
|
||||
Halfedge_around_vertex_const_circulator;
|
||||
|
||||
typedef typename Geometry_traits_2::Point_2 Point_2;
|
||||
typedef typename Geometry_traits_2::Ray_2 Ray_2;
|
||||
typedef typename Geometry_traits_2::Segment_2 Segment_2;
|
||||
typedef typename Geometry_traits_2::Line_2 Line_2;
|
||||
typedef typename Geometry_traits_2::Vector_2 Vector_2;
|
||||
typedef typename Geometry_traits_2::Direction_2 Direction_2;
|
||||
typedef typename Geometry_traits_2::FT Number_type;
|
||||
typedef typename Geometry_traits_2::Object_2 Object_2;
|
||||
|
||||
typedef RegularizationTag Regularization_tag;
|
||||
typedef CGAL::Tag_false Supports_general_polygon_tag;
|
||||
typedef CGAL::Tag_true Supports_simple_polygon_tag;
|
||||
|
||||
Simple_polygon_visibility_2() : p_arr(NULL), geom_traits(NULL) {};
|
||||
|
||||
/*! Constructor given an arrangement and the Regularization tag. */
|
||||
Simple_polygon_visibility_2(const Arrangement_2& arr):
|
||||
p_arr(&arr) {
|
||||
geom_traits = p_arr->geometry_traits();
|
||||
query_pt_is_vertex = false;
|
||||
query_pt_is_on_halfedge = false;
|
||||
}
|
||||
|
||||
|
||||
const std::string name(){return std::string("S_visibility_2");}
|
||||
|
||||
/*! Method to check if the visibility object is attached or not to
|
||||
an arrangement*/
|
||||
bool is_attached() {
|
||||
return (p_arr != NULL);
|
||||
}
|
||||
|
||||
/*! Attaches the visibility object to the 'arr' arrangement */
|
||||
void attach(const Arrangement_2& arr) {
|
||||
p_arr = &arr;
|
||||
geom_traits = p_arr->geometry_traits();
|
||||
query_pt_is_vertex = false;
|
||||
query_pt_is_on_halfedge = false;
|
||||
}
|
||||
|
||||
/*! Detaches the visibility object from the arrangement it is
|
||||
attached to*/
|
||||
void detach() {
|
||||
p_arr = NULL;
|
||||
geom_traits = NULL;
|
||||
vertices.clear();
|
||||
query_pt_is_vertex = false;
|
||||
query_pt_is_on_halfedge = false;
|
||||
p_cdt = boost::shared_ptr<CDT>();
|
||||
}
|
||||
|
||||
/*! Getter method for the input arrangement*/
|
||||
const Arrangement_2& arr() {
|
||||
return *p_arr;
|
||||
}
|
||||
|
||||
/*! Computes the visibility object from the query point 'q' in the face
|
||||
'face' and constructs the output in 'out_arr'*/
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(const Point_2& q,
|
||||
const Face_const_handle face,
|
||||
VARR& out_arr) {
|
||||
|
||||
assert(query_pt_is_vertex == false);
|
||||
assert(query_pt_is_on_halfedge == false);
|
||||
|
||||
// Now retrieve the circulator to first visible vertex from triangulation
|
||||
Ccb_halfedge_const_circulator circ = find_visible_start(face, q);
|
||||
Ccb_halfedge_const_circulator curr = circ;
|
||||
Halfedge_const_handle he;
|
||||
|
||||
do {
|
||||
he = curr;
|
||||
vertices.push_back(he->source()->point());
|
||||
} while(++curr != circ);
|
||||
|
||||
vertices.push_back(vertices[0]);
|
||||
|
||||
visibility_region_impl(q);
|
||||
|
||||
typename std::vector<Point_2> points;
|
||||
while (!s.empty()) {
|
||||
Point_2 curr_point = s.top();
|
||||
points.push_back(curr_point);
|
||||
s.pop();
|
||||
}
|
||||
|
||||
std::reverse(points.begin(), points.end());
|
||||
|
||||
CGAL::Visibility_2::report_while_handling_needles
|
||||
<Simple_polygon_visibility_2>(geom_traits,
|
||||
q,
|
||||
points,
|
||||
out_arr);
|
||||
|
||||
CGAL_precondition(out_arr.number_of_isolated_vertices() == 0);
|
||||
CGAL_precondition(s.size() == 0);
|
||||
conditional_regularize(out_arr, Regularization_tag());
|
||||
vertices.clear();
|
||||
if (out_arr.faces_begin()->is_unbounded()) {
|
||||
return ++out_arr.faces_begin();
|
||||
}
|
||||
else {
|
||||
return out_arr.faces_begin();
|
||||
}
|
||||
}
|
||||
|
||||
/*! Computes the visibility region of the query point 'q' located on the
|
||||
halfedge 'he' and constructs the output in 'out_arr'*/
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(
|
||||
const Point_2& q,
|
||||
const Halfedge_const_handle he,
|
||||
VARR& out_arr )
|
||||
{
|
||||
query_pt_is_vertex = false;
|
||||
query_pt_is_on_halfedge = false;
|
||||
bool query_on_target = false;
|
||||
|
||||
if (q != he->source()->point()) {
|
||||
if (q != he->target()->point()) {
|
||||
vertices.push_back(he->target()->point());
|
||||
query_pt_is_on_halfedge = true;
|
||||
}
|
||||
else {
|
||||
query_pt_is_vertex = true;
|
||||
query_on_target = true;
|
||||
}
|
||||
} else {
|
||||
vertices.push_back( he->target()->point() );
|
||||
query_pt_is_vertex = true;
|
||||
}
|
||||
|
||||
Ccb_halfedge_const_circulator circ = he;
|
||||
circ++;
|
||||
Ccb_halfedge_const_circulator curr = circ;
|
||||
|
||||
do {
|
||||
Halfedge_const_handle he_handle = curr;
|
||||
Point_2 curr_vertex = he_handle->target()->point();
|
||||
vertices.push_back(curr_vertex);
|
||||
} while (++curr != circ);
|
||||
|
||||
if ( query_on_target ) {
|
||||
vertices.push_back( vertices[0] );
|
||||
}
|
||||
|
||||
visibility_region_impl(q);
|
||||
|
||||
typename std::vector<Point_2> points;
|
||||
if (!s.empty()) {
|
||||
Point_2 prev_pt = s.top();
|
||||
if (prev_pt != q) {
|
||||
points.push_back(prev_pt);
|
||||
}
|
||||
else if (query_pt_is_vertex) {
|
||||
points.push_back(prev_pt);
|
||||
}
|
||||
if (!s.empty()) {
|
||||
s.pop();
|
||||
}
|
||||
while(!s.empty()) {
|
||||
Point_2 curr_pt = s.top();
|
||||
if (curr_pt != q) {
|
||||
points.push_back(curr_pt);
|
||||
}
|
||||
else if (query_pt_is_vertex) {
|
||||
points.push_back(curr_pt);
|
||||
}
|
||||
s.pop();
|
||||
}
|
||||
}
|
||||
|
||||
std::reverse(points.begin(), points.end());
|
||||
|
||||
CGAL::Visibility_2::report_while_handling_needles
|
||||
<Simple_polygon_visibility_2>(geom_traits,
|
||||
q,
|
||||
points,
|
||||
out_arr);
|
||||
CGAL_precondition(out_arr.number_of_isolated_vertices() == 0);
|
||||
CGAL_precondition(s.size() == 0);
|
||||
conditional_regularize(out_arr, Regularization_tag());
|
||||
vertices.clear();
|
||||
if (out_arr.faces_begin()->is_unbounded()) {
|
||||
return ++out_arr.faces_begin();
|
||||
}
|
||||
else {
|
||||
return out_arr.faces_begin();
|
||||
}
|
||||
}
|
||||
|
||||
private:
|
||||
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
|
||||
typedef CGAL::Constrained_triangulation_face_base_2<K> Fb;
|
||||
typedef CGAL::Triangulation_data_structure_2<Vb,Fb> TDS;
|
||||
typedef CGAL::No_intersection_tag Itag;
|
||||
typedef CGAL::Constrained_triangulation_2<K, TDS, Itag> CDT;
|
||||
|
||||
private:
|
||||
const Arrangement_2 *p_arr;
|
||||
/*! Boost pointer to the constrained Delaunay triangulation object*/
|
||||
boost::shared_ptr<CDT> p_cdt;
|
||||
/*! Mapping of the vertices of the input to the corresponding circulator
|
||||
needed for finding the first visible vertex in case of face queries*/
|
||||
std::map<Point_2, typename Arrangement_2::Ccb_halfedge_const_circulator>
|
||||
point_itr_map;
|
||||
const Geometry_traits_2 *geom_traits;
|
||||
/*! Stack of visibile points; manipulated when going through the sequence
|
||||
of input vertices; contains the vertices of the visibility region after
|
||||
the run of the algorithm*/
|
||||
std::stack<Point_2> s;
|
||||
/*! Sequence of input vertices*/
|
||||
std::vector<Point_2> vertices;
|
||||
/*! State of visibility region algorithm*/
|
||||
enum {LEFT, RIGHT, SCANA, SCANB, SCANC, SCAND, FINISH} upcase;
|
||||
bool query_pt_is_vertex;
|
||||
bool query_pt_is_on_halfedge;
|
||||
|
||||
/*! Regularize output if flag is set to true*/
|
||||
template <typename VARR>
|
||||
void conditional_regularize(VARR& out_arr, CGAL::Tag_true) {
|
||||
regularize_output(out_arr);
|
||||
}
|
||||
/*! No need to regularize output if flag is set to false*/
|
||||
template <typename VARR>
|
||||
void conditional_regularize(VARR& out_arr, CGAL::Tag_false) {
|
||||
//do nothing
|
||||
}
|
||||
|
||||
/*! Regularizes the output - removes edges that have the same face on both
|
||||
sides */
|
||||
template <typename VARR>
|
||||
void regularize_output(VARR& out_arr) {
|
||||
typename VARR::Edge_iterator e_itr;
|
||||
for (e_itr = out_arr.edges_begin() ;
|
||||
e_itr != out_arr.edges_end() ; e_itr++) {
|
||||
|
||||
typename VARR::Halfedge_handle he = e_itr;
|
||||
typename VARR::Halfedge_handle he_twin = he->twin();
|
||||
if (he->face() == he_twin->face()) {
|
||||
out_arr.remove_edge(he);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/*! Initialized the constrained Delaunay triangulation using the edges of
|
||||
the outer boundary of 'face' */
|
||||
void init_cdt(const Face_const_handle &face) {
|
||||
|
||||
std::vector<std::pair<Point_2,Point_2> > constraints;
|
||||
typename Arrangement_2::Ccb_halfedge_const_circulator circ =
|
||||
face->outer_ccb();
|
||||
typename Arrangement_2::Ccb_halfedge_const_circulator curr = circ;
|
||||
typename Arrangement_2::Halfedge_const_handle he;
|
||||
|
||||
do {
|
||||
he = curr;
|
||||
Point_2 source = he->source()->point();
|
||||
Point_2 target = he->target()->point();
|
||||
point_itr_map.insert(std::make_pair(source, curr));
|
||||
constraints.push_back(std::make_pair(source,target));
|
||||
} while(++curr != circ);
|
||||
|
||||
p_cdt = boost::shared_ptr<CDT>(new CDT(constraints.begin(),constraints.end()));
|
||||
}
|
||||
|
||||
|
||||
/*! Finds a visible vertex from the query point 'q' in 'face'
|
||||
to start the algorithm from*/
|
||||
Ccb_halfedge_const_circulator find_visible_start(Face_const_handle face, const Point_2 &q) {
|
||||
init_cdt(face);
|
||||
typename CDT::Face_handle fh = p_cdt->locate(q);
|
||||
Point_2 start_point = fh->vertex(0)->point();
|
||||
|
||||
// Now retrieve the circulator to first visible vertex from triangulation
|
||||
Ccb_halfedge_const_circulator circ = point_itr_map[start_point];
|
||||
Halfedge_const_handle he_curr = circ;
|
||||
|
||||
Halfedge_around_vertex_const_circulator incident_circ = he_curr->source()->incident_halfedges();
|
||||
Halfedge_around_vertex_const_circulator incident_curr = incident_circ;
|
||||
|
||||
do {
|
||||
Ccb_halfedge_const_circulator curr_inc = incident_curr;
|
||||
Halfedge_const_handle he_curr_inc = curr_inc;
|
||||
|
||||
if (he_curr_inc->face() == face) {
|
||||
Ccb_halfedge_const_circulator incident_next = incident_curr;
|
||||
incident_next++;
|
||||
Halfedge_const_handle he_next_inc = incident_next;
|
||||
|
||||
if (CGAL::Visibility_2::orientation_2(geom_traits,
|
||||
he_curr_inc->source()->point(),
|
||||
he_curr_inc->target()->point(),
|
||||
q) == CGAL::LEFT_TURN
|
||||
|| CGAL::Visibility_2::orientation_2(geom_traits,
|
||||
he_next_inc->source()->point(),
|
||||
he_next_inc->target()->point(),
|
||||
q) == CGAL::LEFT_TURN) {
|
||||
Ccb_halfedge_const_circulator result_circ = incident_next;
|
||||
Halfedge_const_handle he_print = result_circ;
|
||||
return result_circ;
|
||||
}
|
||||
}
|
||||
} while (++incident_curr != incident_circ);
|
||||
}
|
||||
|
||||
/*! Main method of the algorithm - initializes the stack and variables
|
||||
and calles the corresponding methods acc. to the algorithm's state;
|
||||
'q' - query point;
|
||||
'i' - current vertex' index
|
||||
'w' - endpoint of ray shot from query point */
|
||||
void visibility_region_impl(const Point_2& q) {
|
||||
int i = 0;
|
||||
Point_2 w;
|
||||
CGAL::Orientation orient = CGAL::Visibility_2::orientation_2(geom_traits,
|
||||
q,
|
||||
vertices[0],
|
||||
vertices[1]);
|
||||
if ( orient != CGAL::RIGHT_TURN ) {
|
||||
upcase = LEFT;
|
||||
i = 1;
|
||||
w = vertices[1];
|
||||
s.push(vertices[0]);
|
||||
s.push(vertices[1]);
|
||||
}
|
||||
else {
|
||||
upcase = SCANA;
|
||||
i = 1;
|
||||
w = vertices[1];
|
||||
s.push(vertices[0]);
|
||||
}
|
||||
|
||||
Ray_2 ray_origin( q, vertices[0] );
|
||||
do {
|
||||
switch(upcase) {
|
||||
case LEFT:
|
||||
left(i, w, q);
|
||||
break;
|
||||
case RIGHT:
|
||||
right(i, w, q);
|
||||
break;
|
||||
case SCANA:
|
||||
scana(i, w, q);
|
||||
break;
|
||||
case SCANB:
|
||||
scanb(i, w, q);
|
||||
break;
|
||||
case SCANC:
|
||||
scanc(i, w, q);
|
||||
break;
|
||||
case SCAND:
|
||||
scand(i, w, q);
|
||||
break;
|
||||
}
|
||||
if ( upcase == LEFT ) {
|
||||
Point_2 s_t = s.top();
|
||||
s.pop();
|
||||
if ( ( CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, q, vertices[0], s.top() ) == CGAL::RIGHT_TURN ) &&
|
||||
( CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, q, vertices[0],s_t ) == CGAL::LEFT_TURN ) ) {
|
||||
Segment_2 seg( s.top(), s_t );
|
||||
if ( CGAL::Visibility_2::do_intersect_2
|
||||
<Geometry_traits_2, Segment_2, Ray_2>
|
||||
( geom_traits, seg, ray_origin ) ) {
|
||||
Object_2 result = CGAL::Visibility_2::intersect_2
|
||||
<Geometry_traits_2, Segment_2, Ray_2>
|
||||
( geom_traits, seg, ray_origin );
|
||||
const Point_2 * ipoint = CGAL::object_cast<Point_2>(&result);
|
||||
assert( ipoint != NULL );
|
||||
s_t = *ipoint;
|
||||
upcase = SCANB;
|
||||
}
|
||||
}
|
||||
s.push( s_t );
|
||||
}
|
||||
} while(upcase != FINISH);
|
||||
}
|
||||
|
||||
/*! Method that handles the left turns in the vertex algorithm */
|
||||
void left(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
if (i >= vertices.size() - 1) {
|
||||
upcase = FINISH;
|
||||
}
|
||||
else {
|
||||
Point_2 s_t = s.top();
|
||||
s.pop();
|
||||
Point_2 s_t_prev = s.top();
|
||||
s.push( s_t );
|
||||
CGAL::Orientation orient1 = CGAL::Visibility_2::orientation_2
|
||||
<Geometry_traits_2>
|
||||
( geom_traits,
|
||||
query_pt,
|
||||
vertices[i],
|
||||
vertices[i+1] );
|
||||
if ( orient1 != CGAL::RIGHT_TURN ) {
|
||||
// Case L2
|
||||
upcase = LEFT;
|
||||
s.push( vertices[i+1] );
|
||||
w = vertices[i+1];
|
||||
i++;
|
||||
} else {
|
||||
CGAL::Orientation orient2 = CGAL::Visibility_2::orientation_2
|
||||
<Geometry_traits_2>
|
||||
( geom_traits,
|
||||
s_t_prev,
|
||||
vertices[i],
|
||||
vertices[i+1] );
|
||||
if ( orient2 == CGAL::RIGHT_TURN ) {
|
||||
// Case L3
|
||||
upcase = SCANA;
|
||||
w = vertices[i+1];
|
||||
i++;
|
||||
} else {
|
||||
// Case L4
|
||||
upcase = RIGHT;
|
||||
w = vertices[i];
|
||||
i++;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/*! Scans the stack such that all vertices that were pushed before to the
|
||||
stack and are now not visible anymore. */
|
||||
void right(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
Point_2 s_j;
|
||||
Point_2 s_j_prev;
|
||||
Point_2 u;
|
||||
int mode = 0;
|
||||
CGAL::Orientation orient1, orient2;
|
||||
|
||||
s_j_prev = s.top();
|
||||
orient2 = CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, query_pt, s_j_prev, vertices[i] );
|
||||
while ( s.size() > 1 ) {
|
||||
s_j = s_j_prev;
|
||||
orient1 = orient2;
|
||||
s.pop();
|
||||
s_j_prev = s.top();
|
||||
|
||||
orient2 = CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, query_pt, s_j_prev, vertices[i] );
|
||||
if ( orient1 != CGAL::LEFT_TURN && orient2 != CGAL::RIGHT_TURN ) {
|
||||
mode = 1;
|
||||
break;
|
||||
}
|
||||
|
||||
Segment_2 seg2( vertices[i-1], vertices[i] );
|
||||
Segment_2 seg( s_j_prev, s_j );
|
||||
if ( ( vertices[i-1] != s_j )
|
||||
&& ( CGAL::Visibility_2::do_intersect_2
|
||||
<Geometry_traits_2, Segment_2, Segment_2>
|
||||
(geom_traits, seg, seg2) ) ) {
|
||||
Object_2 result = CGAL::Visibility_2::intersect_2
|
||||
<Geometry_traits_2, Segment_2, Segment_2>( geom_traits, seg, seg2 );
|
||||
const Point_2 * ipoint = CGAL::object_cast<Point_2>(&result);
|
||||
assert( ipoint != NULL );
|
||||
u = *ipoint;
|
||||
mode = 2;
|
||||
break;
|
||||
}
|
||||
}
|
||||
|
||||
assert( mode != 0 );
|
||||
if ( mode == 1 ) {
|
||||
orient1 = CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, query_pt, vertices[i], vertices[i+1] );
|
||||
orient2 = CGAL::Visibility_2::orientation_2 <Geometry_traits_2>
|
||||
( geom_traits, vertices[i-1], vertices[i], vertices[i+1] );
|
||||
if ( orient1 == CGAL::RIGHT_TURN ) {
|
||||
// Case R1
|
||||
// Since the next action is RIGHT, we do not compute the intersection
|
||||
// of (s_j,s_j_prev) and the ray (query_pt, vertices[i]),
|
||||
// thus, (s_j,s_j_prev) is not shortcutted, but it is harmless
|
||||
upcase = RIGHT;
|
||||
s.push( s_j );
|
||||
w = vertices[i];
|
||||
i++;
|
||||
} else if ( orient2 == CGAL::RIGHT_TURN ) {
|
||||
// Case R2
|
||||
Ray_2 ray( query_pt, vertices[i] );
|
||||
Segment_2 seg( s_j_prev, s_j );
|
||||
Object_2 result = CGAL::Visibility_2::intersect_2
|
||||
<Geometry_traits_2, Segment_2, Ray_2>
|
||||
( geom_traits, seg, ray );
|
||||
const Point_2 * ipoint = CGAL::object_cast<Point_2>(&result);
|
||||
assert( ipoint != NULL );
|
||||
u = *ipoint;
|
||||
if ( s.top() != u ) {
|
||||
s.push( u );
|
||||
}
|
||||
upcase = LEFT;
|
||||
s.push( vertices[i] );
|
||||
s.push( vertices[i+1] );
|
||||
w = vertices[i+1];
|
||||
i++;
|
||||
} else {
|
||||
// Case R3
|
||||
Ray_2 ray( query_pt, vertices[i] );
|
||||
Segment_2 seg( s_j_prev, s_j );
|
||||
Object_2 result = CGAL::Visibility_2::intersect_2
|
||||
<Geometry_traits_2, Segment_2, Ray_2>
|
||||
( geom_traits, seg, ray );
|
||||
const Point_2 * ipoint = CGAL::object_cast<Point_2>(&result);
|
||||
assert( ipoint != NULL );
|
||||
u = *ipoint;
|
||||
if ( s.top() != u ) {
|
||||
s.push( u );
|
||||
}
|
||||
upcase = SCANC;
|
||||
w = vertices[i];
|
||||
i++;
|
||||
}
|
||||
} else if ( mode == 2 ) {
|
||||
// Case R4
|
||||
upcase = SCAND;
|
||||
w = u;
|
||||
}
|
||||
}
|
||||
|
||||
/*! Scans the vertices starting from index 'i' for the first visible vertex
|
||||
out of the back hidden window */
|
||||
void scana(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
// Scan v_i, v_i+1, ..., v_n for the first edge to intersect (z, s_t)
|
||||
Point_2 u;
|
||||
int k = scan_edges( i, query_pt, s.top(), u, true );
|
||||
|
||||
CGAL::Orientation orient1 = CGAL::Visibility_2::orientation_2
|
||||
<Geometry_traits_2>
|
||||
( geom_traits,
|
||||
query_pt,
|
||||
vertices[k],
|
||||
vertices[k+1] );
|
||||
if ( orient1 == CGAL::RIGHT_TURN ) {
|
||||
bool fwd = CGAL::Visibility_2::collinear_are_ordered_along_line_2
|
||||
<Geometry_traits_2>
|
||||
( geom_traits, query_pt, s.top(), u );
|
||||
if ( !fwd ) {
|
||||
// Case A1
|
||||
upcase = RIGHT;
|
||||
i = k+1;
|
||||
w = u;
|
||||
} else {
|
||||
// Case A2
|
||||
upcase = SCAND;
|
||||
i = k+1;
|
||||
w = u;
|
||||
}
|
||||
} else {
|
||||
// Case A3
|
||||
upcase = LEFT;
|
||||
i = k+1;
|
||||
s.push( u );
|
||||
if ( u != vertices[k+1] ) {
|
||||
s.push( vertices[k+1] );
|
||||
}
|
||||
w = vertices[k+1];
|
||||
}
|
||||
}
|
||||
|
||||
/*! Find the first edge interecting the segment (v_0, s_t) */
|
||||
void scanb(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
if ( i == vertices.size() - 1 ) {
|
||||
upcase = FINISH;
|
||||
return;
|
||||
}
|
||||
Point_2 u;
|
||||
int k = scan_edges( i, s.top(), vertices[0], u, false );
|
||||
if ( (k+1 == vertices.size()-1) && (vertices[0] == u) ) {
|
||||
// Case B1
|
||||
upcase = FINISH;
|
||||
s.push( vertices[0] );
|
||||
} else {
|
||||
// Case B2
|
||||
upcase = RIGHT;
|
||||
i = k+1;
|
||||
w = u;
|
||||
}
|
||||
}
|
||||
|
||||
/*! Finds the exit from a general front hidden window by finding the first
|
||||
vertex to the right of the ray defined by the query_point and w*/
|
||||
void scanc(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
Point_2 u;
|
||||
int k = scan_edges( i, s.top(), w, u, false );
|
||||
upcase = RIGHT;
|
||||
i = k+1;
|
||||
w = u;
|
||||
}
|
||||
|
||||
/*! find the first edge intersecting the given window (s_t, w) */
|
||||
void scand(int& i, Point_2& w, const Point_2& query_pt) {
|
||||
Point_2 u;
|
||||
int k = scan_edges( i, s.top(), w, u, false );
|
||||
upcase = LEFT;
|
||||
i = k+1;
|
||||
s.push( u );
|
||||
if ( u != vertices[k+1] ) {
|
||||
s.push( vertices[k+1] );
|
||||
}
|
||||
w = vertices[k+1];
|
||||
}
|
||||
|
||||
/*! Scan edges v_i,v_{i+1},...,v_n, until find an edge intersecting given ray
|
||||
or given segment. is_ray = true -> ray, false -> segment.
|
||||
The intersection point is returned by u */
|
||||
int scan_edges( int i, const Point_2& ray_begin, const Point_2& ray_end, Point_2& u, bool is_ray ) {
|
||||
CGAL::Orientation old_orient = CGAL::RIGHT_TURN;
|
||||
Ray_2 ray( ray_begin, ray_end );
|
||||
Segment_2 s2( ray_begin, ray_end );
|
||||
int k;
|
||||
Object_2 result;
|
||||
for ( k = i; k+1 < vertices.size(); k++ ) {
|
||||
CGAL::Orientation curr_orient = CGAL::Visibility_2::orientation_2
|
||||
( geom_traits,
|
||||
ray_begin,
|
||||
ray_end,
|
||||
vertices[k+1] );
|
||||
if ( curr_orient != old_orient ) {
|
||||
// Orientation switch, an intersection may occur
|
||||
Segment_2 seg( vertices[k], vertices[k+1] );
|
||||
if ( is_ray ) {
|
||||
if (CGAL::Visibility_2::do_intersect_2
|
||||
<Geometry_traits_2, Segment_2, Ray_2>
|
||||
(geom_traits, seg, ray) ) {
|
||||
result = CGAL::Visibility_2::intersect_2
|
||||
< Geometry_traits_2, Segment_2, Ray_2 >
|
||||
( geom_traits, seg, ray );
|
||||
break;
|
||||
}
|
||||
} else {
|
||||
if (CGAL::Visibility_2::do_intersect_2
|
||||
<Geometry_traits_2, Segment_2, Segment_2>
|
||||
(geom_traits, seg, s2) ) {
|
||||
result = CGAL::Visibility_2::intersect_2
|
||||
< Geometry_traits_2, Segment_2, Segment_2 >
|
||||
( geom_traits, seg, s2 );
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
old_orient = curr_orient;
|
||||
}
|
||||
assert( k+1<vertices.size() );
|
||||
const Point_2 * ipoint = CGAL::object_cast<Point_2>( &result );
|
||||
if ( ipoint ) {
|
||||
u = *ipoint;
|
||||
} else {
|
||||
u = vertices[k+1];
|
||||
}
|
||||
return k;
|
||||
}
|
||||
};
|
||||
|
||||
} // namespace CGAL
|
||||
#endif
|
||||
|
|
@ -0,0 +1,579 @@
|
|||
// Copyright (c) 2013 Technical University Braunschweig (Germany).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s): Michael Hemmer <michael.hemmer@cgal.org>
|
||||
//
|
||||
|
||||
#ifndef CGAL_TRIANGULAR_EXPANSION_VISIBILITY_2__H
|
||||
#define CGAL_TRIANGULAR_EXPANSION_VISIBILITY_2__H
|
||||
|
||||
#include <CGAL/Arrangement_2.h>
|
||||
#include <boost/shared_ptr.hpp>
|
||||
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
|
||||
|
||||
namespace CGAL {
|
||||
|
||||
template<class Arrangement_2_ , class RegularizationTag = CGAL::Tag_true >
|
||||
class Triangular_expansion_visibility_2 {
|
||||
typedef typename Arrangement_2_::Geometry_traits_2 Geometry_traits_2;
|
||||
typedef typename Geometry_traits_2::Kernel K;
|
||||
public:
|
||||
// Currently only consider with same type for both
|
||||
typedef Arrangement_2_ Arrangement_2;
|
||||
typedef typename Arrangement_2::Traits_2 Traits_2;
|
||||
typedef typename Arrangement_2::Halfedge_const_handle Halfedge_const_handle;
|
||||
typedef typename Arrangement_2::Halfedge_handle Halfedge_handle;
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator
|
||||
Ccb_halfedge_const_circulator;
|
||||
typedef typename Arrangement_2::Face_const_handle Face_const_handle;
|
||||
typedef typename Arrangement_2::Face_handle Face_handle;
|
||||
typedef typename Arrangement_2::Vertex_const_handle Vertex_const_handle;
|
||||
typedef typename Arrangement_2::Vertex_handle Vertex_handle;
|
||||
|
||||
typedef typename K::Point_2 Point_2;
|
||||
typedef typename Geometry_traits_2::Ray_2 Ray_2;
|
||||
typedef typename Geometry_traits_2::Segment_2 Segment_2;
|
||||
typedef typename Geometry_traits_2::Line_2 Line_2;
|
||||
typedef typename Geometry_traits_2::Vector_2 Vector_2;
|
||||
typedef typename Geometry_traits_2::Direction_2 Direction_2;
|
||||
typedef typename Geometry_traits_2::FT Number_type;
|
||||
typedef typename Geometry_traits_2::Object_2 Object_2;
|
||||
|
||||
// TODO
|
||||
typedef RegularizationTag Regularization_tag;
|
||||
|
||||
typedef CGAL::Tag_true Supports_general_polygon_tag;
|
||||
typedef CGAL::Tag_true Supports_simple_polygon_tag;
|
||||
|
||||
private:
|
||||
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
|
||||
typedef CGAL::Constrained_triangulation_face_base_2<K> Fb;
|
||||
typedef CGAL::Triangulation_data_structure_2<Vb,Fb> TDS;
|
||||
typedef CGAL::No_intersection_tag Itag;
|
||||
typedef CGAL::Constrained_Delaunay_triangulation_2<K, TDS, Itag> CDT;
|
||||
|
||||
public:
|
||||
const std::string name(){return std::string("T_visibility_2");}
|
||||
private:
|
||||
const Arrangement_2* p_arr;
|
||||
boost::shared_ptr<CDT> p_cdt;
|
||||
std::vector<Segment_2> needles;
|
||||
|
||||
public:
|
||||
Triangular_expansion_visibility_2() : p_arr(NULL){}
|
||||
|
||||
/*! Constructor given an arrangement and the Regularization tag. */
|
||||
Triangular_expansion_visibility_2 (const Arrangement_2& arr)
|
||||
: p_arr(&arr){
|
||||
init_cdt();
|
||||
}
|
||||
|
||||
bool is_attached() {
|
||||
//std::cout << "is_attached" << std::endl;
|
||||
return (p_arr != NULL);
|
||||
}
|
||||
|
||||
void attach(const Arrangement_2& arr) {
|
||||
// todo observe changes in arr;
|
||||
if(p_arr != &arr){
|
||||
p_arr = &arr;
|
||||
init_cdt();
|
||||
}
|
||||
//std::cout << "attach done" << std::endl;
|
||||
}
|
||||
|
||||
void detach() {
|
||||
//std::cout << "detach" << std::endl;
|
||||
p_arr = NULL;
|
||||
p_cdt = boost::shared_ptr<CDT>();
|
||||
}
|
||||
|
||||
const Arrangement_2& arr() {
|
||||
return *p_arr;
|
||||
}
|
||||
|
||||
typename CDT::Edge get_edge(typename CDT::Face_handle fh, int i){
|
||||
return std::make_pair(fh,i);
|
||||
}
|
||||
|
||||
Point_2 ray_seg_intersection(
|
||||
const Point_2& q, const Point_2& b, // the ray
|
||||
const Point_2& s, const Point_2& t) // the segment
|
||||
{
|
||||
Ray_2 ray(q,b);
|
||||
Segment_2 seg(s,t);
|
||||
assert(typename K::Do_intersect_2()(ray,seg));
|
||||
CGAL::Object obj = typename K::Intersect_2()(ray,seg);
|
||||
Point_2 result = object_cast<Point_2>(obj);
|
||||
return result;
|
||||
}
|
||||
|
||||
void collect_needle(
|
||||
const Point_2& q,
|
||||
const typename CDT::Vertex_handle vh,
|
||||
const typename CDT::Face_handle fh,
|
||||
int index){
|
||||
|
||||
|
||||
// the expanded edge should not be constrained
|
||||
assert(!p_cdt->is_constrained(get_edge(fh,index)));
|
||||
assert(!p_cdt->is_infinite(fh));
|
||||
// go into the new face
|
||||
const typename CDT::Face_handle nfh(fh->neighbor(index));
|
||||
assert(!p_cdt->is_infinite(nfh));
|
||||
|
||||
// get indices of neighbors
|
||||
int nindex = nfh->index(fh); // index of new vertex and old face
|
||||
int rindex = p_cdt->ccw(nindex); // index of face behind right edge
|
||||
int lindex = p_cdt-> cw(nindex); // index of face behind left edge
|
||||
|
||||
// get vertices seen from entering edge
|
||||
const typename CDT::Vertex_handle nvh(nfh->vertex(nindex));
|
||||
const typename CDT::Vertex_handle rvh(nfh->vertex(p_cdt->cw (nindex)));
|
||||
const typename CDT::Vertex_handle lvh(nfh->vertex(p_cdt->ccw(nindex)));
|
||||
assert(!p_cdt->is_infinite(nvh));
|
||||
assert(!p_cdt->is_infinite(lvh));
|
||||
assert(!p_cdt->is_infinite(rvh));
|
||||
|
||||
// get edges seen from entering edge
|
||||
typename CDT::Edge re = get_edge(nfh,p_cdt->ccw(nindex));
|
||||
typename CDT::Edge le = get_edge(nfh,p_cdt-> cw(nindex));
|
||||
|
||||
// do orientation computation once for new vertex
|
||||
typename K::Orientation_2 orientation =
|
||||
p_cdt->geom_traits().orientation_2_object();
|
||||
CGAL::Orientation orient = orientation(q,vh->point(),nvh->point());
|
||||
|
||||
|
||||
//std::cout << "\n collect_needle" <<std::endl;
|
||||
//std::cout << "q "<< q << std::endl ;
|
||||
//std::cout << "vh->point() "<< vh->point() << std::endl;
|
||||
//std::cout << "lvh->point() "<< lvh->point() << std::endl ;
|
||||
//std::cout << "nvh->point() "<< nvh->point() << std::endl ;
|
||||
//std::cout << "rvh->point() "<< rvh->point() << std::endl<< std::endl;
|
||||
|
||||
|
||||
|
||||
switch ( orient ) {
|
||||
case CGAL::COUNTERCLOCKWISE:
|
||||
// looking on to the right edge
|
||||
if(p_cdt->is_constrained(re)){
|
||||
if(vh!=rvh){
|
||||
Point_2 p = ray_seg_intersection(q,vh->point(),nvh->point(),rvh->point());
|
||||
//std::cout << vh->point() <<" -1- "<< p <<std::endl;
|
||||
needles.push_back(Segment_2(vh->point(),p));
|
||||
}
|
||||
}else{
|
||||
collect_needle(q,vh,nfh,rindex);
|
||||
}
|
||||
break;
|
||||
case CGAL::CLOCKWISE:
|
||||
// looking on to the left edge
|
||||
if(p_cdt->is_constrained(le)){
|
||||
if(vh!=lvh){
|
||||
Point_2 p = ray_seg_intersection(q,vh->point(),nvh->point(),lvh->point());
|
||||
//std::cout << vh->point() <<" -2- "<< p <<std::endl;
|
||||
needles.push_back(Segment_2(vh->point(),p));
|
||||
}
|
||||
}else{
|
||||
collect_needle(q,vh,nfh,lindex);
|
||||
}
|
||||
break;
|
||||
default:
|
||||
assert(orient == CGAL::COLLINEAR);
|
||||
// looking on nvh, so it must be reported
|
||||
// if it wasn't already (triangles rotate around vh)
|
||||
if(vh != nvh){
|
||||
//std::cout << vh->point() <<" -3- "<< nvh->point() <<std::endl;
|
||||
needles.push_back(Segment_2(vh->point(),nvh->point()));
|
||||
}
|
||||
// but we may also contiue looking along the vertex
|
||||
if(!p_cdt->is_constrained(re)){
|
||||
collect_needle(q,nvh,nfh,rindex);
|
||||
}
|
||||
if(!p_cdt->is_constrained(le)){
|
||||
collect_needle(q,nvh,nfh,lindex);
|
||||
}
|
||||
break;
|
||||
}
|
||||
}
|
||||
template<class OIT>
|
||||
OIT expand_edge(
|
||||
const Point_2& q,
|
||||
const Point_2& left,
|
||||
const Point_2& right,
|
||||
typename CDT::Face_handle fh,
|
||||
int index,
|
||||
OIT oit){
|
||||
|
||||
// the expanded edge should not be constrained
|
||||
assert(!p_cdt->is_constrained(get_edge(fh,index)));
|
||||
assert(!p_cdt->is_infinite(fh));
|
||||
|
||||
// go into the new face
|
||||
const typename CDT::Face_handle nfh(fh->neighbor(index));
|
||||
assert(!p_cdt->is_infinite(nfh));
|
||||
|
||||
// get indices of neighbors
|
||||
int nindex = nfh->index(fh); // index of new vertex and old face
|
||||
int rindex = p_cdt->ccw(nindex); // index of face behind right edge
|
||||
int lindex = p_cdt-> cw(nindex); // index of face behind left edge
|
||||
|
||||
// get vertices seen from entering edge
|
||||
const typename CDT::Vertex_handle nvh(nfh->vertex(nindex));
|
||||
const typename CDT::Vertex_handle rvh(nfh->vertex(p_cdt->cw (nindex)));
|
||||
const typename CDT::Vertex_handle lvh(nfh->vertex(p_cdt->ccw(nindex)));
|
||||
assert(!p_cdt->is_infinite(nvh));
|
||||
assert(!p_cdt->is_infinite(lvh));
|
||||
assert(!p_cdt->is_infinite(rvh));
|
||||
|
||||
// get edges seen from entering edge
|
||||
typename CDT::Edge re = get_edge(nfh,p_cdt->ccw(nindex));
|
||||
typename CDT::Edge le = get_edge(nfh,p_cdt-> cw(nindex));
|
||||
|
||||
// do orientation computation once for new vertex
|
||||
typename K::Orientation_2 orientation =
|
||||
p_cdt->geom_traits().orientation_2_object();
|
||||
CGAL::Orientation ro = orientation(q,right,nvh->point());
|
||||
CGAL::Orientation lo = orientation(q,left ,nvh->point());
|
||||
|
||||
assert(typename K::Orientation_2()(q,left ,lvh->point()) != CGAL::CLOCKWISE);
|
||||
assert(typename K::Orientation_2()(q,right,rvh->point()) != CGAL::COUNTERCLOCKWISE);
|
||||
|
||||
//std::cout << (ro == CGAL::COUNTERCLOCKWISE) << " " << (lo == CGAL::CLOCKWISE) << std::endl;
|
||||
|
||||
|
||||
|
||||
//right edge is seen if new vertex is counter clockwise of right boarder
|
||||
if(ro == CGAL::COUNTERCLOCKWISE){
|
||||
if(p_cdt->is_constrained(re)){
|
||||
// the edge is constrained
|
||||
// report intersection with right boarder ray
|
||||
// if it is not already the right vertex (already reported)
|
||||
if(right != rvh->point()){
|
||||
*oit++ = ray_seg_intersection(q,right,nvh->point(),rvh->point());
|
||||
}
|
||||
|
||||
// then report intersection with left boarder if it exists
|
||||
if(lo == CGAL::COUNTERCLOCKWISE){
|
||||
*oit++ = ray_seg_intersection(q,left,nvh->point(),rvh->point());
|
||||
}
|
||||
}else{
|
||||
// the edge is not a constrained
|
||||
if(lo == CGAL::COUNTERCLOCKWISE){
|
||||
// no split needed and return
|
||||
//std::cout<< "h1"<< std::endl;
|
||||
oit = expand_edge(q,left,right,nfh,rindex,oit);
|
||||
//std::cout<< "h1 done"<< std::endl;
|
||||
return oit;
|
||||
}else{
|
||||
// spliting at new vertex
|
||||
//std::cout<< "h2"<< std::endl;
|
||||
*oit++ = expand_edge(q,nvh->point(),right,nfh,rindex,oit);
|
||||
//std::cout<< "h2 done"<< std::endl;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
//std::cout << "q "<< q << std::endl ;
|
||||
//std::cout << "lvh->point() "<< lvh->point() << std::endl;
|
||||
//std::cout << "left "<< left << std::endl ;
|
||||
//std::cout << "nvh->point() "<< nvh->point() << std::endl ;
|
||||
//std::cout << "right "<< right << std::endl ;
|
||||
//std::cout << "rvh->point() "<< rvh->point() << std::endl<< std::endl;
|
||||
|
||||
|
||||
// determin whether new vertex needs to be reported
|
||||
if(ro != CGAL::CLOCKWISE && lo != CGAL::COUNTERCLOCKWISE){
|
||||
*oit++ = nvh->point();
|
||||
}
|
||||
if(!Regularization_tag::value){
|
||||
assert(!(ro == CGAL::COLLINEAR && lo == CGAL::COLLINEAR));
|
||||
// we have to check whether a needle starts here.
|
||||
if(p_cdt->is_constrained(le) && !p_cdt->is_constrained(re) && ro == CGAL::COLLINEAR)
|
||||
collect_needle(q,nvh,nfh,rindex);
|
||||
if(p_cdt->is_constrained(re) && !p_cdt->is_constrained(le) && lo == CGAL::COLLINEAR)
|
||||
collect_needle(q,nvh,nfh,lindex);
|
||||
}
|
||||
|
||||
//left edge is seen if new vertex is clockwise of left boarder
|
||||
if(lo == CGAL::CLOCKWISE){
|
||||
if(p_cdt->is_constrained(le)){
|
||||
// the edge is constrained
|
||||
// report interesection with right boarder if exists
|
||||
if(ro == CGAL::CLOCKWISE){
|
||||
*oit++ = ray_seg_intersection(q,right,nvh->point(),lvh->point());
|
||||
}
|
||||
// then report intersection with left boarder ray
|
||||
// if it is not already the left vertex (already reported)
|
||||
if(left != lvh->point()){
|
||||
*oit++ = ray_seg_intersection(q,left,nvh->point(),lvh->point());
|
||||
}
|
||||
return oit;
|
||||
}else{
|
||||
// the edge is not a constrained
|
||||
if(ro == CGAL::CLOCKWISE){
|
||||
// no split needed and return
|
||||
//std::cout<< "h3"<< std::endl;
|
||||
oit = expand_edge(q,left,right,nfh,lindex,oit);
|
||||
//std::cout<< "h3 done"<< std::endl;
|
||||
return oit;
|
||||
}else{
|
||||
// spliting at new vertex
|
||||
//std::cout<< "h4"<< std::endl;
|
||||
oit = expand_edge(q,left,nvh->point(),nfh,lindex,oit);
|
||||
//std::cout<< "h4 done"<< std::endl;
|
||||
return oit;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
return oit;
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(const Point_2& q,
|
||||
const Face_const_handle face,
|
||||
VARR& out_arr
|
||||
){
|
||||
//std::cout << "query in face interior" << std::endl;
|
||||
|
||||
out_arr.clear();
|
||||
needles.clear();
|
||||
assert(!face->is_unbounded());
|
||||
|
||||
|
||||
std::vector<Point_2> raw_output;
|
||||
typename CDT::Face_handle fh = p_cdt->locate(q);
|
||||
|
||||
raw_output.push_back(fh->vertex(1)->point());
|
||||
if(!p_cdt->is_constrained(get_edge(fh,0))){
|
||||
//std::cout<< "edge 0 is not constrained" << std::endl;
|
||||
expand_edge(
|
||||
q,
|
||||
fh->vertex(2)->point(),
|
||||
fh->vertex(1)->point(),
|
||||
fh,0,std::back_inserter(raw_output));
|
||||
}
|
||||
|
||||
raw_output.push_back(fh->vertex(2)->point());
|
||||
if(!p_cdt->is_constrained(get_edge(fh,1))){
|
||||
//std::cout << "edge 1 is not constrained" << std::endl;
|
||||
expand_edge(
|
||||
q,
|
||||
fh->vertex(0)->point(),
|
||||
fh->vertex(2)->point(),
|
||||
fh,1,std::back_inserter(raw_output));
|
||||
}
|
||||
|
||||
raw_output.push_back(fh->vertex(0)->point());
|
||||
if(!p_cdt->is_constrained(get_edge(fh,2))){
|
||||
//std::cout << "edge 2 is not constrained" << std::endl;
|
||||
expand_edge(
|
||||
q,
|
||||
fh->vertex(1)->point(),
|
||||
fh->vertex(0)->point(),
|
||||
fh,2,std::back_inserter(raw_output));
|
||||
}
|
||||
|
||||
|
||||
return output(raw_output,out_arr);
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
compute_visibility(const Point_2& q,
|
||||
const Halfedge_const_handle he,
|
||||
VARR& out_arr) {
|
||||
//std::cout << "visibility_region he" << std::endl;
|
||||
|
||||
assert(!he->face()->is_unbounded());
|
||||
out_arr.clear();
|
||||
needles.clear();
|
||||
|
||||
std::vector<Point_2> raw_output;
|
||||
typename CDT::Locate_type location;
|
||||
int index;
|
||||
typename CDT::Face_handle fh = p_cdt->locate(q,location,index);
|
||||
assert(location == CDT::EDGE || location == CDT::VERTEX);
|
||||
// the following code tries to figure out which triangle one should start in.
|
||||
|
||||
|
||||
if(location == CDT::EDGE){
|
||||
//std::cout << "query on edge" << std::endl;
|
||||
// this is the easy part, there are only two possible faces
|
||||
// index indicates the edge = vertex on the other side of the edge
|
||||
// the next vertex in cw order should be the target of given edge
|
||||
if(fh->vertex(p_cdt->cw(index))->point() != he->target()->point()){
|
||||
//std::cout << "need to swap face" << std::endl;
|
||||
// take face on the other side if this is not the case
|
||||
typename CDT::Face_handle nfh = fh->neighbor(index);
|
||||
index = nfh->index(fh);
|
||||
fh = nfh;
|
||||
}
|
||||
assert(fh->vertex(p_cdt->cw(index))->point() == he->target()->point());
|
||||
assert(!p_cdt->is_infinite(fh->vertex(index)));
|
||||
|
||||
|
||||
// output the edge the query lies on
|
||||
raw_output.push_back(he->source()->point());
|
||||
raw_output.push_back(he->target()->point());
|
||||
|
||||
if(!p_cdt->is_constrained(get_edge(fh,p_cdt->ccw(index)))){
|
||||
expand_edge(
|
||||
q,
|
||||
fh->vertex(index)->point(), //left
|
||||
he->target()->point() , //right
|
||||
fh,
|
||||
p_cdt->ccw(index),
|
||||
std::back_inserter(raw_output));
|
||||
}
|
||||
raw_output.push_back(fh->vertex(index)->point());
|
||||
|
||||
if(!p_cdt->is_constrained(get_edge(fh,p_cdt->cw(index)))){
|
||||
expand_edge(
|
||||
q,
|
||||
he->source()->point() , //left
|
||||
fh->vertex(index)->point(), //right
|
||||
fh,
|
||||
p_cdt->cw(index),
|
||||
std::back_inserter(raw_output));
|
||||
}
|
||||
}
|
||||
|
||||
if(location == CDT::VERTEX){
|
||||
//std::cout << "query on vertex" << std::endl;
|
||||
|
||||
//bool query_point_on_vertex_is_not_working_yet = false;
|
||||
//assert(query_point_on_vertex_is_not_working_yet);
|
||||
|
||||
assert(q == he->target()->point());
|
||||
assert(fh->vertex(index)->point() == he->target()->point());
|
||||
|
||||
// push points that are seen anyway
|
||||
// raw_output.push_back(he->source()->point()); inserted last
|
||||
raw_output.push_back(he->target()->point());
|
||||
raw_output.push_back(he->next()->target()->point());
|
||||
|
||||
// now start in the triangle that contains he->next()
|
||||
while(
|
||||
p_cdt->is_infinite(fh->vertex(p_cdt->ccw(index))) ||
|
||||
he->next()->target()->point() != fh->vertex(p_cdt->ccw(index))->point()){
|
||||
typename CDT::Face_handle nfh = fh->neighbor(p_cdt->ccw(index));
|
||||
int nindex = nfh->index(fh);
|
||||
index = p_cdt->ccw(nindex);
|
||||
fh = nfh;
|
||||
assert(he->target()->point() == fh->vertex(index)->point());
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
assert(he->next()->source()->point() == fh->vertex(index)->point());
|
||||
assert(he->next()->target()->point() == fh->vertex(p_cdt->ccw(index))->point());
|
||||
assert(!p_cdt->is_infinite(fh));
|
||||
assert(p_cdt->is_constrained(get_edge(fh,p_cdt->cw(index))));
|
||||
|
||||
while(he->source()->point() != fh->vertex(p_cdt->ccw(index))->point()){
|
||||
|
||||
if(!p_cdt->is_constrained(get_edge(fh,index))){
|
||||
expand_edge(
|
||||
q,
|
||||
fh->vertex(p_cdt-> cw(index))->point(), //left
|
||||
fh->vertex(p_cdt->ccw(index))->point(), //right
|
||||
fh,
|
||||
index,
|
||||
std::back_inserter(raw_output));
|
||||
}
|
||||
// push left end point of edge into output
|
||||
raw_output.push_back(fh->vertex(p_cdt-> cw(index))->point());
|
||||
|
||||
// take the next triangle around q in ccw order
|
||||
typename CDT::Face_handle nfh = fh->neighbor(p_cdt->ccw(index));
|
||||
int nindex = nfh->index(fh);
|
||||
index = p_cdt->ccw(nindex);
|
||||
fh = nfh;
|
||||
assert(fh->vertex(index)->point() == he->target()->point());
|
||||
}
|
||||
}
|
||||
return output(raw_output,out_arr);
|
||||
}
|
||||
|
||||
template <typename VARR>
|
||||
typename VARR::Face_handle
|
||||
output(std::vector<Point_2>& raw_output, VARR& out_arr){
|
||||
|
||||
if(needles.size()>0){
|
||||
std::vector<Segment_2> segments(needles.begin(),needles.end());
|
||||
for(unsigned int i = 0; i <raw_output.size();i++){
|
||||
// //std::cout << raw_output[i] << " -- "
|
||||
// << raw_output[(i+1)%raw_output.size()] << std::endl;
|
||||
segments.push_back(Segment_2(raw_output[i],raw_output[(i+1)%raw_output.size()]));
|
||||
}
|
||||
|
||||
CGAL::insert_non_intersecting_curves(out_arr,segments.begin(),segments.end());
|
||||
//CGAL::insert(out_arr,segments.begin(),segments.end());
|
||||
}else{
|
||||
typename VARR::Vertex_handle v_last, v_first;
|
||||
v_last = v_first =
|
||||
out_arr.insert_in_face_interior(raw_output[0], out_arr.unbounded_face());
|
||||
|
||||
for(unsigned int i = 0; i <raw_output.size()-1;i++){
|
||||
// //std::cout << raw_output[i] << " -- "
|
||||
// << raw_output[(i+1)%raw_output.size()] << std::endl;
|
||||
if(raw_output[i]<raw_output[(i+1)]){
|
||||
v_last = out_arr.insert_from_left_vertex (Segment_2(raw_output[i], raw_output[i+1]), v_last)->target();
|
||||
}else{
|
||||
v_last = out_arr.insert_from_right_vertex(Segment_2(raw_output[i], raw_output[i+1]), v_last)->target();
|
||||
}
|
||||
}
|
||||
out_arr.insert_at_vertices(Segment_2(raw_output.front(),raw_output.back()),v_last,v_first);
|
||||
}
|
||||
|
||||
assert(out_arr.number_of_faces()== 2);
|
||||
|
||||
if(out_arr.faces_begin()->is_unbounded())
|
||||
return ++out_arr.faces_begin();
|
||||
else
|
||||
return out_arr.faces_begin();
|
||||
}
|
||||
|
||||
void init_cdt(){
|
||||
//std::cout<< "==============" <<std::endl;
|
||||
//std::cout<< "Input Polygon:" <<std::endl;
|
||||
//todo, avoid copy by using modified iterator
|
||||
std::vector<std::pair<Point_2,Point_2> > constraints;
|
||||
for(typename Arrangement_2::Edge_const_iterator eit = p_arr->edges_begin();
|
||||
eit != p_arr->edges_end(); eit++){
|
||||
Point_2 source = eit->source()->point();
|
||||
Point_2 target = eit->target()->point();
|
||||
//std::cout << source << " -- " << target << std::endl;
|
||||
constraints.push_back(std::make_pair(source,target));
|
||||
}
|
||||
//std::cout << "init_cdt new CDT" << std::endl;
|
||||
p_cdt = boost::shared_ptr<CDT>(new CDT(constraints.begin(),constraints.end()));
|
||||
//std::cout << "init_cdt done" << std::endl;
|
||||
//std::cout << std::endl;
|
||||
}
|
||||
};
|
||||
|
||||
} // namespace CGAL
|
||||
|
||||
#endif //CGAL_TRIANGULAR_EXPANSION_VISIBILITY_2__H
|
||||
|
|
@ -0,0 +1,371 @@
|
|||
// Copyright (c) 2013 Technical University Braunschweig (Germany).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org).
|
||||
// You can redistribute it and/or modify it under the terms of the GNU
|
||||
// General Public License as published by the Free Software Foundation,
|
||||
// either version 3 of the License, or (at your option) any later version.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
//
|
||||
// Author(s): Francisc Bungiu <fbungiu@gmail.com>
|
||||
// Michael Hemmer <michael.hemmer@cgal.org>
|
||||
|
||||
#ifndef CGAL_VISIBILITY_UTILS_H
|
||||
#define CGAL_VISIBILITY_UTILS_H
|
||||
|
||||
#include <vector>
|
||||
#include <CGAL/tags.h>
|
||||
#include <CGAL/enum.h>
|
||||
|
||||
namespace CGAL {
|
||||
namespace Visibility_2 {
|
||||
|
||||
template <class Arrangement_2>
|
||||
int count_edges_in_face(typename Arrangement_2::Face_const_handle fch) {
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator
|
||||
Ccb_halfedge_const_circulator;
|
||||
|
||||
Ccb_halfedge_const_circulator circ = fch->outer_ccb();
|
||||
Ccb_halfedge_const_circulator curr = circ;
|
||||
|
||||
int edge_cnt(0);
|
||||
do {
|
||||
edge_cnt++;
|
||||
} while (++curr != circ);
|
||||
return edge_cnt;
|
||||
}
|
||||
|
||||
template <class Edge_const_iterator>
|
||||
void print_edge(Edge_const_iterator eit) {
|
||||
std::cout << "[" << eit->curve() << "]" << std::endl;
|
||||
}
|
||||
template <class Face_const_handle, class Ccb_halfedge_const_circulator>
|
||||
void print_simple_face(Face_const_handle fh) {
|
||||
Ccb_halfedge_const_circulator cir = fh->outer_ccb();
|
||||
Ccb_halfedge_const_circulator curr = cir;
|
||||
do {
|
||||
std::cout << "[" << curr->curve() << "]" << std::endl;
|
||||
} while (++ curr != cir);
|
||||
}
|
||||
|
||||
template <class Arrangement_2>
|
||||
void print_arrangement(const Arrangement_2& arr) {
|
||||
typedef typename Arrangement_2::Edge_const_iterator Edge_const_iterator;
|
||||
Edge_const_iterator eit;
|
||||
std::cout << arr.number_of_edges() << " edges:" << std::endl;
|
||||
for (eit = arr.edges_begin(); eit != arr.edges_end(); ++eit)
|
||||
print_edge(eit);
|
||||
}
|
||||
|
||||
template <class Arrangement_2>
|
||||
void print_arrangement_by_face(const Arrangement_2& arr) {
|
||||
typedef typename Arrangement_2::Face_const_iterator Face_const_iterator;
|
||||
typedef typename Arrangement_2::Ccb_halfedge_const_circulator
|
||||
Ccb_halfedge_const_circulator;
|
||||
Face_const_iterator fit;
|
||||
for (fit = arr.faces_begin() ; fit != arr.faces_end() ; fit++) {
|
||||
if (!fit->is_unbounded()) {
|
||||
std::cout << "FACE\n";
|
||||
print_simple_face<Face_const_iterator, Ccb_halfedge_const_circulator>(fit);
|
||||
std::cout << "END FACE\n";
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2>
|
||||
Orientation orientation_2(const Geometry_traits_2 *geom_traits,
|
||||
const typename Geometry_traits_2::Point_2& p,
|
||||
const typename Geometry_traits_2::Point_2& q,
|
||||
const typename Geometry_traits_2::Point_2& r) {
|
||||
|
||||
typename Geometry_traits_2::Orientation_2 orient =
|
||||
geom_traits->orientation_2_object();
|
||||
return orient(p, q, r);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2>
|
||||
bool less_distance_to_point_2(const Geometry_traits_2 *geom_traits,
|
||||
const typename Geometry_traits_2::Point_2& p,
|
||||
const typename Geometry_traits_2::Point_2& q,
|
||||
const typename Geometry_traits_2::Point_2& r) {
|
||||
|
||||
typename Geometry_traits_2::Less_distance_to_point_2 less_dist =
|
||||
geom_traits->less_distance_to_point_2_object();
|
||||
return less_dist(p, q, r);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2>
|
||||
bool collinear(const Geometry_traits_2 *geom_traits,
|
||||
const typename Geometry_traits_2::Point_2& p,
|
||||
const typename Geometry_traits_2::Point_2& q,
|
||||
const typename Geometry_traits_2::Point_2& r) {
|
||||
|
||||
typename Geometry_traits_2::Collinear_2 collinear_fnct =
|
||||
geom_traits->collinear_2_object();
|
||||
return collinear_fnct(p, q, r);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2, class _Curve_first, class _Curve_second >
|
||||
typename Geometry_traits_2::Object_2 intersect_2(const Geometry_traits_2 *geom_traits,
|
||||
const _Curve_first& s1,
|
||||
const _Curve_second& s2) {
|
||||
|
||||
typedef typename Geometry_traits_2::Kernel Kernel;
|
||||
const Kernel *kernel = static_cast<const Kernel*> (geom_traits);
|
||||
typename Kernel::Intersect_2 intersect_fnct = kernel->intersect_2_object();
|
||||
return intersect_fnct(s1, s2);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2>
|
||||
CGAL::Comparison_result compare_xy_2(const Geometry_traits_2 *geom_traits,
|
||||
const typename Geometry_traits_2::Point_2 &p,
|
||||
const typename Geometry_traits_2::Point_2 &q) {
|
||||
|
||||
typename Geometry_traits_2::Compare_xy_2 cmp =
|
||||
geom_traits->compare_xy_2_object();
|
||||
return cmp(p, q);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2, class Type1, class Type2>
|
||||
typename Geometry_traits_2::FT compute_squared_distance_2(
|
||||
const Geometry_traits_2 *geom_traits,
|
||||
const Type1& p,
|
||||
const Type2& seg) {
|
||||
|
||||
typename Geometry_traits_2::Compute_squared_distance_2 compute_dist =
|
||||
geom_traits->compute_squared_distance_2_object();
|
||||
return compute_dist(p, seg);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2, class Type1, class Type2>
|
||||
bool do_intersect_2(const Geometry_traits_2 *geom_traits,
|
||||
const Type1& c1,
|
||||
const Type2& c2) {
|
||||
|
||||
typename Geometry_traits_2::Do_intersect_2 intersect =
|
||||
geom_traits->do_intersect_2_object();
|
||||
return intersect(c1, c2);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2>
|
||||
bool collinear_are_ordered_along_line_2(const Geometry_traits_2 *geom_traits,
|
||||
const typename Geometry_traits_2::Point_2 &p,
|
||||
const typename Geometry_traits_2::Point_2 &q,
|
||||
const typename Geometry_traits_2::Point_2 &r) {
|
||||
|
||||
typename Geometry_traits_2::Collinear_are_ordered_along_line_2 coll =
|
||||
geom_traits->collinear_are_ordered_along_line_2_object();
|
||||
return coll(p, q, r);
|
||||
}
|
||||
|
||||
template <class Geometry_traits_2, class Type1, class Type2>
|
||||
typename Geometry_traits_2::Point_2 construct_projected_point_2(
|
||||
const Geometry_traits_2 *geom_traits,
|
||||
const Type1& s,
|
||||
const Type2& p) {
|
||||
|
||||
typedef typename Geometry_traits_2::Point_2 Point_2;
|
||||
typedef typename Geometry_traits_2::FT Number_type;
|
||||
typename Geometry_traits_2::Construct_projected_point_2 construct_proj =
|
||||
geom_traits->construct_projected_point_2_object();
|
||||
Point_2 proj_pt = construct_proj(s.supporting_line(), p);
|
||||
if (s.has_on(proj_pt)) {
|
||||
return proj_pt;
|
||||
}
|
||||
else {
|
||||
Number_type d_to_src = compute_squared_distance_2
|
||||
<Geometry_traits_2, Point_2, Point_2>(geom_traits, proj_pt, s.source());
|
||||
Number_type d_to_trg = compute_squared_distance_2
|
||||
<Geometry_traits_2, Point_2, Point_2>(geom_traits, proj_pt, s.target());
|
||||
if (d_to_src < d_to_trg) {
|
||||
return s.source();
|
||||
}
|
||||
else {
|
||||
return s.target();
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
//construct an arrangement of visibility region from a vector of circular ordered vertices with respect to the query point
|
||||
template <class Visibility_2, class Visibility_arrangement_2>
|
||||
void report_while_handling_needles(
|
||||
const typename Visibility_2::Arrangement_2::Geometry_traits_2 *geom_traits,
|
||||
const typename Visibility_2::Arrangement_2::Point_2& q,
|
||||
std::vector<typename Visibility_2::Arrangement_2::Point_2>& points,
|
||||
Visibility_arrangement_2& arr_out) {
|
||||
|
||||
typedef typename Visibility_2::Arrangement_2 Arrangement_2;
|
||||
typedef typename Arrangement_2::Point_2 Point_2;
|
||||
typedef typename Arrangement_2::Geometry_traits_2 Geometry_traits_2;
|
||||
typedef typename Arrangement_2::Halfedge_handle Halfedge_handle;
|
||||
typedef typename Arrangement_2::Vertex_handle Vertex_handle;
|
||||
typedef typename Geometry_traits_2::Segment_2 Segment_2;
|
||||
typedef typename Geometry_traits_2::Direction_2 Direction_2;
|
||||
|
||||
typename std::vector<Segment_2>::size_type i = 0;
|
||||
|
||||
if (points.front() == points.back()) {
|
||||
points.pop_back();
|
||||
}
|
||||
|
||||
while (collinear(geom_traits,
|
||||
q,
|
||||
points[i],
|
||||
points.back())) {
|
||||
|
||||
points.push_back(points[i]);
|
||||
i++;
|
||||
}
|
||||
|
||||
points.push_back(points[i]);
|
||||
|
||||
typename Visibility_arrangement_2::Halfedge_handle he_handle;
|
||||
typename Visibility_arrangement_2::Vertex_handle v_trg; //the handle of vertex where the next segment is inserted
|
||||
typename Visibility_arrangement_2::Vertex_handle v_fst; //the handle of vertex inserted first
|
||||
typename Visibility_arrangement_2::Vertex_handle v_needle_end; //the handle of vertex of the end of a needle
|
||||
|
||||
v_trg = v_fst = arr_out.insert_in_face_interior(points[i], arr_out.unbounded_face());
|
||||
//find a point that is right after a needle
|
||||
while (i+1 < points.size()) {
|
||||
if ( collinear(geom_traits,
|
||||
points[i],
|
||||
points[i+1],
|
||||
q)) {
|
||||
typename Visibility_arrangement_2::Vertex_handle v_needle_begin = v_trg;
|
||||
std::vector<Point_2> forward_needle; //vertices of the needle that are not between q and v_needle_begin; their direction is leaving q;
|
||||
std::vector<Point_2> backward_needle; //vertices of the needle that are not between q and v_needle_begin; their direction is towards q;
|
||||
std::vector<Point_2> part_in_q_side; //vertices of the needle that are between q and v_needle_begin
|
||||
part_in_q_side.push_back(points[i]);
|
||||
forward_needle.push_back((points[i]));
|
||||
|
||||
bool same_side_of_q = (compare_xy_2(geom_traits, points[i], q)==compare_xy_2(geom_traits, points[i], points[i+1]));
|
||||
if (same_side_of_q)
|
||||
part_in_q_side.push_back(points[i+1]);
|
||||
else
|
||||
forward_needle.push_back(points[i+1]);
|
||||
i++;
|
||||
while (i+1< points.size() && orientation_2(geom_traits,
|
||||
points[i],
|
||||
points[i+1],
|
||||
q ) == CGAL::COLLINEAR) {
|
||||
if (same_side_of_q) {
|
||||
part_in_q_side.push_back(points[i+1]);
|
||||
}
|
||||
else {
|
||||
if (compare_xy_2(geom_traits, part_in_q_side.front(), q)
|
||||
== compare_xy_2(geom_traits, part_in_q_side.front(), points[i+1])) {
|
||||
same_side_of_q = true;
|
||||
part_in_q_side.push_back(points[i+1]);
|
||||
}
|
||||
else {
|
||||
if (less_distance_to_point_2(geom_traits, q, points[i], points[i+1]))
|
||||
forward_needle.push_back(points[i+1]);
|
||||
else
|
||||
backward_needle.push_back(points[i+1]);
|
||||
}
|
||||
}
|
||||
i++;
|
||||
}
|
||||
//obtain the end point of a needle
|
||||
Point_2 end_of_needle;
|
||||
if (same_side_of_q)
|
||||
end_of_needle = part_in_q_side.back();
|
||||
else {
|
||||
if (backward_needle.empty()) {
|
||||
end_of_needle = forward_needle.back();
|
||||
}
|
||||
else {
|
||||
end_of_needle = backward_needle.back();
|
||||
}
|
||||
}
|
||||
|
||||
std::reverse(backward_needle.begin(), backward_needle.end());
|
||||
std::vector<Point_2> merged_needle;
|
||||
|
||||
// merge the forward_needle and backward_needle
|
||||
unsigned int itr_fst = 0, itr_snd = 0;
|
||||
while (itr_fst < forward_needle.size() &&
|
||||
itr_snd < backward_needle.size()) {
|
||||
|
||||
if (less_distance_to_point_2(geom_traits,
|
||||
q,
|
||||
forward_needle[itr_fst],
|
||||
backward_needle[itr_snd])) {
|
||||
|
||||
merged_needle.push_back(forward_needle[itr_fst]);
|
||||
itr_fst++;
|
||||
}
|
||||
else {
|
||||
merged_needle.push_back(backward_needle[itr_snd]);
|
||||
itr_snd++;
|
||||
}
|
||||
}
|
||||
while (itr_fst < forward_needle.size()) {
|
||||
merged_needle.push_back(forward_needle[itr_fst]);
|
||||
itr_fst++;
|
||||
}
|
||||
while (itr_snd < backward_needle.size()) {
|
||||
merged_needle.push_back(backward_needle[itr_snd]);
|
||||
itr_snd++;
|
||||
}
|
||||
|
||||
for (unsigned int p = 0 ; p+1 < merged_needle.size() ; p++) {
|
||||
if (CGAL::Visibility_2::compare_xy_2<Geometry_traits_2>(geom_traits, merged_needle[p], merged_needle[p+1]) == CGAL::SMALLER) {
|
||||
he_handle = arr_out.insert_from_left_vertex(Segment_2(merged_needle[p], merged_needle[p+1]), v_trg);
|
||||
}
|
||||
else {
|
||||
he_handle = arr_out.insert_from_right_vertex(Segment_2(merged_needle[p], merged_needle[p+1]), v_trg);
|
||||
}
|
||||
v_trg = he_handle->target();
|
||||
if (merged_needle[p+1] == end_of_needle) {
|
||||
v_needle_end = v_trg;
|
||||
}
|
||||
}
|
||||
if (same_side_of_q) {
|
||||
//insert the part of needle between q and v_needle_begin
|
||||
v_trg = v_needle_begin;
|
||||
for (unsigned int p = 0 ; p+1 < part_in_q_side.size() ; p++) {
|
||||
if (CGAL::Visibility_2::compare_xy_2<Geometry_traits_2>(geom_traits, part_in_q_side[p], part_in_q_side[p+1]) == CGAL::SMALLER) {
|
||||
he_handle = arr_out.insert_from_left_vertex(Segment_2(part_in_q_side[p], part_in_q_side[p+1]), v_trg);
|
||||
}
|
||||
else {
|
||||
he_handle = arr_out.insert_from_right_vertex(Segment_2(part_in_q_side[p], part_in_q_side[p+1]), v_trg);
|
||||
}
|
||||
v_trg = he_handle->target();
|
||||
}
|
||||
}
|
||||
else
|
||||
v_trg = v_needle_end;
|
||||
}
|
||||
else {
|
||||
if (CGAL::Visibility_2::compare_xy_2<Geometry_traits_2>(geom_traits, v_trg->point(), points[i+1]) == CGAL::SMALLER) {
|
||||
he_handle = arr_out.insert_from_left_vertex(Segment_2(points[i], points[i+1]), v_trg);
|
||||
}
|
||||
else {
|
||||
he_handle = arr_out.insert_from_right_vertex(Segment_2(points[i], points[i+1]), v_trg);
|
||||
}
|
||||
v_trg = he_handle->target();
|
||||
i++;
|
||||
}
|
||||
if (i+2 == points.size()) {
|
||||
//close the boundary
|
||||
v_trg = he_handle->target();
|
||||
arr_out.insert_at_vertices(Segment_2(points[points.size()-2], points[points.size()-1]), v_trg, v_fst);
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
} // end namespace Visibility_2
|
||||
} // end namespace CGAL
|
||||
|
||||
#endif
|
||||
Loading…
Reference in New Issue