segment_delaunay_graph_hierarchy_2_impl.h

很多二维 三维几何计算算法 C++ 类库

  // find the first conflict

  // first look if there are intersections...
  Vertex_circulator vc = hierarchy[0]->incident_vertices(vertices[0]);
  Vertex_circulator vc_start = vc;
  do {
    Vertex_handle vv(vc);
    if ( is_infinite(vv) ) {
      vc++;
      continue;
    }

    Arrangement_type at_res = this->arrangement_type(t, vv);
    if ( vv->is_segment() ) {
      if ( at_res == AT2::DISJOINT || at_res == AT2::TOUCH_1 ||
	   at_res == AT2::TOUCH_2 || at_res == AT2::TOUCH_11 ||
	   at_res == AT2::TOUCH_12 || at_res == AT2::TOUCH_21 ||
	   at_res == AT2::TOUCH_22 ) {
	// do nothing
      } else if ( at_res == AT2::IDENTICAL ) {
	// merge info of identical items
	merge_info(vv, ss);
	return vv;
      } else if ( at_res == AT2::CROSSING ) {
	return insert_intersecting_segment_with_tag(ss, t, vv, level,
						    itag, stag);
      } else if ( at_res == AT2::TOUCH_11_INTERIOR_1 ) {
	Vertex_handle vp = second_endpoint_of_segment(vv);

	Storage_site_2 sss = split_storage_site(ss, vp->storage_site(), true);
	// merge the info of the first (common) subsegment
	merge_info(vv, sss);

	return insert_segment_on_point(ss, vp, level, 1);
      } else if ( at_res == AT2::TOUCH_12_INTERIOR_1 ) {
	Vertex_handle vp = first_endpoint_of_segment(vv);
	return insert_segment_on_point(ss, vp, level, 0);
      } else {
	// this should never be reached; the only possible values for
	// at_res are DISJOINT, CROSSING, TOUCH_11_INTERIOR_1
	// and TOUCH_12_INTERIOR_1
	CGAL_assertion( false );
      }
    } else {
      CGAL_assertion( vv->is_point() );
      if ( at_res == AT2::INTERIOR ) {
	Storage_site_2 svv = vv->storage_site();
	if ( svv.is_input() ) {
	  return insert_segment_on_point(ss, vv, level, 2);
	} else {
	  // MK::ERROR:: not ready yet
	  CGAL_assertion( false );
	}
      }
    }
    ++vc;
  } while ( vc != vc_start );

  // first look for conflict with vertex
  Face_circulator fc_start = hierarchy[0]->incident_faces(vertices[0]);
  Face_circulator fc = fc_start;
  Face_handle start_f;
  Sign s;

  std::map<Face_handle,Sign> sign_map;

  do {
    Face_handle f(fc);

    s = incircle(f, t);

    sign_map[f] = s;

    if ( s == NEGATIVE ) {
      start_f = f;
      break;
    }
    ++fc;
  } while ( fc != fc_start );

  // segments must have a conflict with at least one vertex
  CGAL_assertion( s == NEGATIVE );

  // we are in conflict with a Voronoi vertex; start from that and 
  // find the entire conflict region and then repair the diagram
  List l;
  Face_map fm;

  Triple<bool, Vertex_handle, Arrangement_type>
    vcross(false, Vertex_handle(), AT2::DISJOINT);

  hierarchy[0]->initialize_conflict_region(start_f, l);
  hierarchy[0]->expand_conflict_region(start_f, t, ss, l, fm,
				       sign_map, vcross);

  CGAL_assertion( vcross.third == AT2::DISJOINT ||
		  vcross.third == AT2::CROSSING ||
		  vcross.third == AT2::INTERIOR );

  // the following condition becomes true only if intersecting
  // segments are found
  if ( vcross.first ) {
    if ( t.is_segment() ) {
      if ( vcross.third == AT2::CROSSING ) {
	Intersections_tag itag;
	return insert_intersecting_segment_with_tag(ss, t, vcross.second,
						    level, itag, stag);
      } else if ( vcross.third == AT2::INTERIOR ) {
	return insert_segment_on_point(ss, vcross.second, level, 2);
      } else {
	// this should never be reached; the only possible values for
	// vcross.third are CROSSING, INTERIOR and DISJOINT
	CGAL_assertion( false );
      }
    }
  }

  // no intersecting segment has been found; we insert the segment as
  // usual...
  Vertex_handle v = hierarchy[0]->create_vertex(ss);

  hierarchy[0]->retriangulate_conflict_region(v, l, fm);

  insert_segment_in_upper_levels(t, ss, v, vertices, level, stag);

  return v;
}


template<class Gt, class ST, class STag, class DS, class LTag>
void
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::
insert_segment_in_upper_levels(const Site_2& t, const Storage_site_2& ss,
			       Vertex_handle vbelow,
			       const Vertex_handle* vertices,
			       int level, Tag_true /* stag */)
{
  CGAL_precondition( vertices != NULL );
  CGAL_precondition( vbelow != Vertex_handle() );

  // insert at all upper levels
  Vertex_handle previous = vbelow;
  Vertex_handle vertex = vbelow;

  int k = 1;
  while ( k <= level ) {
    if ( hierarchy[k]->number_of_vertices() == 2 ) {
      Vertex_handle v0(hierarchy[k]->finite_vertices_begin());
      Vertex_handle v1(++(hierarchy[k]->finite_vertices_begin()));
      CGAL_precondition( v0 != Vertex_handle() &&
			 v1 != Vertex_handle() );
      vertex = hierarchy[k]->insert_third(ss, v0, v1);
    } else {
      vertex = hierarchy[k]->insert_segment_interior(t, ss, vertices[k]);
    }

    CGAL_assertion( vertex != Vertex_handle() );

    vertex->set_down(previous); // link with level above
    previous->set_up(vertex);
    previous = vertex;
    k++;
  }
}

//--------------------------------------------------------------------
// insertion of a segment that goes through a point
//--------------------------------------------------------------------

template<class Gt, class ST, class STag, class DS, class LTag>
typename
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::Vertex_handle
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::
insert_segment_on_point(const Storage_site_2& ss,
			const Vertex_handle& v,
			int level, int which)
{  
  // inserts the segment represented by ss in the case where this
  // segment goes through a point which has already been inserted and
  // corresponds to the vertex handle v
  CGAL_precondition( ss.is_segment() );
  CGAL_precondition( v->is_point() );

  Storage_site_2 ssv = v->storage_site();
  Site_2 sv = ssv.site();

  Vertex_handle vnear[sdg_hierarchy_2__maxlevel];

  nearest_neighbor(sv, vnear, false);

  int v_level = find_level(v);

  if ( v_level < level ) {
    Vertex_handle vbelow = vnear[v_level];

    insert_point(sv, ssv, v_level + 1, level, vbelow, vnear);
  }

  // merge the info of the (common) splitting endpoint
  merge_info(v, ss);

  if ( which == 2 ) {
    Storage_site_2 ss1 = this->split_storage_site(ss, ssv, true);
    Storage_site_2 ss2 = this->split_storage_site(ss, ssv, false);

    insert_segment_interior(ss1.site(), ss1, vnear, level);
    return insert_segment_interior(ss2.site(), ss2, vnear, level);
  }

  Storage_site_2 ss1 = this->split_storage_site(ss, ssv, which == 0);
  return insert_segment_interior(ss1.site(), ss1, vnear, level);
}

//--------------------------------------------------------------------
// insertion of an intersecting segment
//--------------------------------------------------------------------
template<class Gt, class ST, class STag, class DS, class LTag>
typename
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::Vertex_handle
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::
insert_intersecting_segment_with_tag(const Storage_site_2& ss,
				     const Site_2& t, Vertex_handle v,
				     int level,
				     Tag_true itag, Tag_false /* stag */)
{
  CGAL_precondition( t.is_segment() && v->is_segment() );

  const Storage_site_2& ssitev = v->storage_site();
  Site_2 sitev = ssitev.site();

  if ( same_segments(t, sitev) ) {
    merge_info(v, ss);
    return v;
  }

  Vertex_triple vt = hierarchy[0]->insert_point_on_segment(ss, t, v, itag);

  Vertex_handle verticesx[sdg_hierarchy_2__maxlevel];

  Vertex_handle vsx = vt.first;
  verticesx[0] = vsx;

  bool compute_new_level = true;
  int new_level = compute_new_level ? random_level() : level;
  if ( new_level > 0 ) {
    Storage_site_2 ssx = vsx->storage_site();
    Site_2 sx = ssx.site();
    insert_point(sx, ssx, 1, new_level, vsx, verticesx);
  }

  Storage_site_2 ss3 =
    this->st_.construct_storage_site_2_object()(ss, ssitev, true);
  Storage_site_2 ss4 =
    this->st_.construct_storage_site_2_object()(ss, ssitev, false);

  Site_2 s3 = ss3.site();
  Site_2 s4 = ss4.site();

  insert_segment_interior(s3, ss3, verticesx, level);
  insert_segment_interior(s4, ss4, verticesx, level);
  return vsx;
}

template<class Gt, class ST, class STag, class DS, class LTag>
typename
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::Vertex_handle
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::
insert_intersecting_segment_with_tag(const Storage_site_2& ss,
				     const Site_2& t, Vertex_handle v,
				     int level,
				     Tag_true itag, Tag_true /* stag */)
{
  CGAL_precondition( t.is_segment() && v->is_segment() );

  CGAL_expensive_precondition( arrangement_type(t, v->site()) );

  const Storage_site_2& ssitev = v->storage_site();
  Site_2 sitev = ssitev.site();

  if ( same_segments(t, sitev) ) {
    // MK::ERROR: I may need to insert it to levels higher than its
    // previous level...
    merge_info(v, ss);
    return v;
  }

  Vertex_handle verticesx[sdg_hierarchy_2__maxlevel];

  int levelv = find_level(v);

  Vertex_handle vcross = v;
  Vertex_handle v1_old, v2_old, vsx_old;

  int k = 0;
  while ( k <= levelv ) {
    Vertex_handle vcross_up = vcross->up();

    Vertex_triple vt =
      hierarchy[k]->insert_point_on_segment(ss, t, vcross, itag);

    // now I need to update the sites for vertices v1 and v2
    Vertex_handle vsx = vt.first;
    Vertex_handle v1 = vt.second;
    Vertex_handle v2 = vt.third;

    CGAL_assertion( v1->is_segment() && v2->is_segment() );

    if ( k > 0 ) {
      CGAL_assertion( (same_segments(v1->site(), v1_old->site()) &&
		       same_segments(v2->site(), v2_old->site())) ||
		      (same_segments(v1->site(), v2_old->site()) &&
		       same_segments(v2->site(), v1_old->site()))
		      );
      if ( same_segments(v1->site(), v1_old->site()) ) {
	v1->set_down(v1_old);
	v2->set_down(v2_old);
	v1_old->set_up(v1);
	v2_old->set_up(v2);
      } else {
	v1->set_down(v2_old);
	v2->set_down(v1_old);
	v1_old->set_up(v2);
	v2_old->set_up(v1);
      }
      vsx_old->set_up(vsx);
      vsx->set_down(vsx_old);
    }

    v1_old = v1;
    v2_old = v2;
    vsx_old = vsx;

    verticesx[k] = vsx;

    vcross = vcross_up;
    k++;
  }

  if ( levelv < level ) {
    Storage_site_2 ssx = verticesx[0]->storage_site();
    Site_2 sx = ssx.site();

    insert_point(sx, ssx, levelv + 1, level, verticesx[levelv], verticesx);
  }

  Storage_site_2 ss3 =
    this->st_.construct_storage_site_2_object()(ss, ssitev, true);

  Storage_site_2 ss4 =
    this->st_.construct_storage_site_2_object()(ss, ssitev, false);

  Site_2 s3 = ss3.site();
  Site_2 s4 = ss4.site();

  insert_segment_interior(s3, ss3, verticesx, level);
  insert_segment_interior(s4, ss4, verticesx, level);
  return verticesx[0];
}


//====================================================================
//====================================================================
//                   METHODS FOR REMOVAL
//====================================================================
//====================================================================

template<class Gt, class ST, class STag, class DS, class LTag>
bool
Segment_Delaunay_graph_hierarchy_2<Gt,ST,STag,DS,LTag>::
remove(const Vertex_handle& v)
{
  CGAL_precondition( !is_infinite(v) );
  const Storage_site_2& ss = v->storage_site();

  if ( !ss.is_input() ) { return false; }

  Point_handle h1, h2;
  bool is_point = ss.is_point();
  if ( is_point ) {
    h1 = ss.point();
  } else {
    CGAL_assertion( ss.is_segment() );   
    h1 = ss.source_of_supporting_site();
    h2 = ss.target_of_supporting_site();