snc_point_locator.h

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

  virtual void transform(const Aff_transformation_3& t) {
    candidate_provider->transform(t);
  }

  virtual bool update( Unique_hash_map<Vertex_handle, bool>& V, 
                       Unique_hash_map<Halfedge_handle, bool>& E, 
                       Unique_hash_map<Halffacet_handle, bool>& F) {
    CGAL_NEF_TIMER(ct_t.start());
    CGAL_assertion( initialized);
    bool updated = candidate_provider->update( V, E, F);
    CGAL_NEF_TIMER(ct_t.stop());
    return updated;
  }

  virtual ~SNC_point_locator_by_spatial_subdivision() {
    CGAL_warning( initialized); // required?
    delete candidate_provider;
  }

  virtual Object_handle shoot(const Ray_3& ray, int mask=255) const {
    CGAL_NEF_TIMER(rs_t.start());
    CGAL_assertion( initialized);
    _CGAL_NEF_TRACEN( "shooting: "<<ray);
    Object_handle result;
    Vertex_handle v;
    Halfedge_handle e;
    Halffacet_handle f;
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
    Halffacet_triangle_handle t;
#endif
#ifdef CGAL_NEF3_FACET_WITH_BOX
    Partial_facet pf;
#endif
    bool hit = false;
    Point_3 eor; // 'end of ray', the latest ray's hit point
    Objects_along_ray objects = candidate_provider->objects_along_ray(ray);
    Objects_along_ray_iterator objects_iterator = objects.begin();
    while( !hit && objects_iterator != objects.end()) {
      Object_list candidates = *objects_iterator;
      Object_list_iterator o;
      CGAL_for_each( o, candidates) {
        if( CGAL::assign( v, *o) && ((mask&1) != 0)) {
          _CGAL_NEF_TRACEN("trying vertex on "<<v->point());
          if( ray.source() != v->point() && ray.has_on(v->point())) {
            _CGAL_NEF_TRACEN("the ray intersects the vertex");
            _CGAL_NEF_TRACEN("prev. intersection? "<<hit);
            if( hit) _CGAL_NEF_TRACEN("prev. intersection on "<<eor);
            if( hit && !Segment_3( ray.source(), eor).has_on(v->point()))
              continue;
            eor = v->point();
            result = Object_handle(v);
            hit = true;
            _CGAL_NEF_TRACEN("the vertex becomes the new hit object");
          }
        }
        else if( CGAL::assign( e, *o) && ((mask&2) != 0)) {
          Point_3 q;
          _CGAL_NEF_TRACEN("trying edge on "<< Segment_3(e->source()->point(),e->twin()->source()->point()));
          if( is.does_intersect_internally( ray, Segment_3(e->source()->point(),
                                                           e->twin()->source()->point()), q)) {
            _CGAL_NEF_TRACEN("ray intersects edge on "<<q);
            _CGAL_NEF_TRACEN("prev. intersection? "<<hit);
            if( hit) _CGAL_NEF_TRACEN("prev. intersection on "<<eor);
            if( hit && !has_smaller_distance_to_point( ray.source(), q, eor))
              continue;
            _CGAL_NEF_TRACEN("is the intersection point on the current cell? "<<
                    candidate_provider->is_point_on_cell( q, objects_iterator));
            if( !candidate_provider->is_point_on_cell( q, objects_iterator))
                continue;
            eor = q; 
            result = Object_handle(e);
            hit = true;
            _CGAL_NEF_TRACEN("the edge becomes the new hit object");
          }
        }
        else if( CGAL::assign( f, *o) && ((mask&4) != 0)) {
          Point_3 q;
          _CGAL_NEF_TRACEN("trying facet with on plane "<<f->plane()<<
                  " with point on "<<f->plane().point());
          if( is.does_intersect_internally( ray, f, q) ) {
            _CGAL_NEF_TRACEN("ray intersects facet on "<<q);
            _CGAL_NEF_TRACEN("prev. intersection? "<<hit);
            if( hit) _CGAL_NEF_TRACEN("prev. intersection on "<<eor);
            if( hit && !has_smaller_distance_to_point( ray.source(), q, eor))
	      continue;
            _CGAL_NEF_TRACEN("is the intersection point on the current cell? "<<
                    candidate_provider->is_point_on_cell( q, objects_iterator));
            if( !candidate_provider->is_point_on_cell( q, objects_iterator))
	      continue;
            eor = q;
            result = Object_handle(f);
            hit = true; 
            _CGAL_NEF_TRACEN("the facet becomes the new hit object");
          }
        }
#ifdef CGAL_NEF3_FACET_WITH_BOX
	else if( CGAL::assign(pf, *o) && ((mask&4) != 0)) {
	  CGAL_NEF_TRACEN("new ray shooting");
          Point_3 q;
          if( is.does_intersect_internally( ray, pf, q) ) {
            if( hit && !has_smaller_distance_to_point( ray.source(), q, eor))
	      continue;
            if( !candidate_provider->is_point_on_cell( q, objects_iterator))
	      continue;
            eor = q;
            result = Object_handle(pf.f);
            hit = true; 
          }
	}
#endif
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
        else if( CGAL::assign( t, *o) && ((mask&8) != 0)) {
          Point_3 q;
          Triangle_3 tr = t.get_triangle();
          _CGAL_NEF_TRACEN("trying triangle "<<tr);
          if( is.does_intersect( ray, tr, q) && normalized(q) != ray.source()) {
            _CGAL_NEF_TRACEN("ray intersect triangle on "<< normalized(q));
            _CGAL_NEF_TRACEN("prev. intersection? "<<hit);
            if( hit) _CGAL_NEF_TRACEN("prev. intersection on "<< normalized(eor));
            if( hit && !has_smaller_distance_to_point( ray.source(), q, eor))
              continue;
            _CGAL_NEF_TRACEN("is the intersection point on the boundary of the facet? "<<
                    is.does_contain_on_boundary( t, q));
            if( is.does_contain_on_boundary( t, q))
              continue;
            _CGAL_NEF_TRACEN("is the intersection point on the current cell? "<<
                    candidate_provider->is_point_on_cell(q,objects_iterator));
            if( !candidate_provider->is_point_on_cell( q, objects_iterator))
                continue;
            eor = q;
            result = Object_handle(Halffacet_handle(t));
            hit = true; 
            _CGAL_NEF_TRACEN("the facet becomes the new hit object");
          }
        }
#endif
        else if((mask&15) == 15)
          CGAL_assertion_msg( 0, "wrong handle");
      }
      if(!hit)
        ++objects_iterator;
    }
    CGAL_NEF_TIMER(rs_t.stop());
    return result;
  }

  virtual Object_handle locate( const Point_3& p) const {
    if(Infi_box::extended_kernel()) {
    CGAL_NEF_TIMER(pl_t.start());
    CGAL_assertion( initialized);
    _CGAL_NEF_TRACEN( "locate "<<p);
    Object_handle result;
    Vertex_handle v;
    Halfedge_handle e;
    Halffacet_handle f;
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
    Halffacet_triangle_handle t;
#endif
#ifdef CGAL_NEF3_FACET_WITH_BOX
    Partial_facet pf;
#endif
    Object_list candidates = candidate_provider->objects_around_point(p);
    Object_list_iterator o = candidates.begin();
    bool found = false;
    while( !found && o != candidates.end()) {
      if( CGAL::assign( v, *o)) {
        if ( p == v->point()) {
          _CGAL_NEF_TRACEN("found on vertex "<<v->point())          
          result = Object_handle(v);
          found = true;
        }
      }
      else if( CGAL::assign( e, *o)) {
        if ( is.does_contain_internally(Segment_3(e->source()->point(),e->twin()->source()->point()), p) ) {
          _CGAL_NEF_TRACEN("found on edge "<<Segment_3(e->source()->point(),e->twin()->source()->point()));
          result = Object_handle(e);
          found = true;
        }
      }
      else if( CGAL::assign( f, *o)) {
        if (is.does_contain_internally( f, p) ) {
          _CGAL_NEF_TRACEN("found on facet...");
          result = Object_handle(f);
          found = true;
        }
      }
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
      else if( CGAL::assign( t, *o)) {
        Triangle_3 tr = t.get_triangle();
        if( tr.has_on(p)) {
          _CGAL_NEF_TRACEN("found on triangle "<<tr);
	  if(is.does_contain_on_boundary( t, p)) {
	    _CGAL_NEF_TRACEN("but located on the facet's boundary");
	    continue;
	  }
          result = Object_handle(Halffacet_handle(t));
	  found = true;
        }
      }
#endif
#ifdef CGAL_NEF3_FACET_WITH_BOX 
      else if( CGAL::assign(pf, *o)) {
	CGAL_NEF_TRACEN("new locate ");
        if ( is.does_contain_internally( pf, p) ) {
          _CGAL_NEF_TRACEN("found on partial facet...");
          result = Object_handle(pf.f);
          found = true;
        }
      }
#endif
      o++;
    }
    if( !found) {
      _CGAL_NEF_TRACEN("point not found in 2-skeleton");
      _CGAL_NEF_TRACEN("shooting ray to determine the volume");
      Ray_3 r( p, Vector_3( -1, 0, 0));
      result = Object_handle(determine_volume(r));
    }    CGAL_NEF_TIMER(pl_t.start());
    CGAL_NEF_TIMER(pl_t.stop());
    return result;
  } else {   // standard kernel
    CGAL_assertion( initialized);
    _CGAL_NEF_TRACEN( "locate "<<p);
    typename SNC_structure::FT min_distance;
    typename SNC_structure::FT tmp_distance;
    Object_handle result;
    Vertex_handle v;
    Halfedge_handle e;
    Halffacet_handle f;
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
    Halffacet_triangle_handle t;
#endif
#ifdef CGAL_NEF3_FACET_WITH_BOX
    Partial_facet pf;
#endif
    Object_list candidates = candidate_provider->objects_around_point(p);
    Object_list_iterator o = candidates.begin();

    if(candidates.empty())
      return Base(*this).volumes_begin();

    CGAL::assign(v,*o);
    CGAL_assertion(CGAL::assign(v,*o));
    if(p==v->point())
      return Object_handle(v);

    min_distance = CGAL::squared_distance(v->point(),p);
    result = Object_handle(v);
    ++o;
    while(o!=candidates.end() && CGAL::assign(v,*o)) {
      if ( p == v->point()) {
        _CGAL_NEF_TRACEN("found on vertex "<<v->point());
        return Object_handle(v);
      }
      tmp_distance = CGAL::squared_distance(v->point(),p);
      if(tmp_distance < min_distance) {
	result = Object_handle(v);
        min_distance = tmp_distance;
      }
      ++o;
    }     

    CGAL::assign(v, result);
    Segment_3 s(p,v->point());
    Point_3 ip;

    /*
    // TODO: das geht effizienter
    Object_list_iterator of(o);
    while(of != candidates.end() && assign(e, *of)) ++of;
    
    typename SNC_structure::SHalfedge_iterator sei;
    for(sei=v->shalfedges_begin(); sei!=v->shalfedges_end(); ++sei){      
      if(sei->is_twin()) continue;
      Halffacet_handle fout = sei->facet();
      if(fout->is_twin()) fout = fout->twin();
      Object_list_iterator ofc(of); 
      for(;ofc!=candidates.end();++ofc) {
	if(CGAL::assign(f,*ofc)) {
	  if(f == fout->twin())
	    std::cerr << "shit" << std::endl;
	  if(f == fout) {
	    Object_list_iterator oe(ofc);
	    --ofc;
	    candidates.erase(oe);
	  }
	}
      }
    }
    */
    for(;o!=candidates.end();++o) {
      if( CGAL::assign( e, *o)) {
	Segment_3 ss(e->source()->point(),e->twin()->source()->point());
	CGAL_NEF_TRACEN("test edge " << e->source()->point() << "->" << e->twin()->source()->point());
	  if(is.does_contain_internally(ss, p) ) {
	_CGAL_NEF_TRACEN("found on edge "<< ss);
          return Object_handle(e);
        }
	if(is.does_intersect_internally(s, ss, ip)) {
	  s = Segment_3(p, normalized(ip));
	  result = Object_handle(e);
        }

      } else 
      if( CGAL::assign( f, *o)) {
	CGAL_NEF_TRACEN("test facet " << f->plane());
        if (is.does_contain_internally(f,p) ) {
          _CGAL_NEF_TRACEN("found on facet...");
          return Object_handle(f);
        }
        if( is.does_intersect_internally(s,f,ip)) {
          s = Segment_3(p, normalized(ip));
	  result = Object_handle(f);
        }
#ifdef CGAL_NEF3_FACET_WITH_BOX
      } else if( CGAL::assign(pf, *o)) {
	CGAL_NEF_TRACEN("new locate ");
	if ( is.does_contain_internally( pf, p) ) {
	  _CGAL_NEF_TRACEN("found on partial facet...");
	  return Object_handle(pf.f);
	}
        if( is.does_intersect_internally(s,pf,ip)) {	
          s = Segment_3(p, normalized(ip));
	  result = Object_handle(pf.f);
        }
#endif
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
      } else if( CGAL::assign(t, *o)) {
        CGAL_NEF_TRACEN("test triangle of facet " << t->plane());
        Triangle_3 tr = t.get_triangle();
        CGAL_NEF_TRACEN("trying triangle "<<tr);	
	if(tr.has_on(p)) {
	  f = t;
	  return Object_handle(f);
        }
        if( is.does_intersect( s, tr, ip)) {
	  CGAL_assertion(!is.does_contain_on_boundary( t, ip));
          s = Segment_3(p, normalized(ip));
	  result = Object_handle(t);
        }
#endif
      } 
      else CGAL_assertion_msg(false, "wrong handle type");
    }

    if( CGAL::assign( v, result)) {
      _CGAL_NEF_TRACEN("vertex hit, obtaining volume...");
      SM_point_locator L(&*v);
      Object_handle so = L.locate(s.source()-s.target());
      SFace_handle sf;
      if(CGAL::assign(sf,so))
        return sf->volume();
      std::cerr << "Abbruch " << std::endl;
      return Object_handle();
      CGAL_assertion_msg(false, "wrong handle type");
/*
      SHalfedge_handle se;
      CGAL_assertion(CGAL::assign(se,so));
      CGAL_NEF_TRACEN("intersect segment " << s << " with edges");
      for(;ox!=candidates.end();++ox) {
	if(!CGAL::assign(e,*ox)) continue;
	CGAL_NEF_TRACEN("test edge " << e->source()->point() << "->" << e->twin()->source()->point());
	if(is.does_intersect_internally(s,Segment_3(e->source()->point(),e->twin()->source()->point()),ip)) {
	  s = Segment_3(p, normalized(ip));
	  result = Object_handle(e);
        }
      }
      CGAL_assertion(CGAL::assign(e,result));
      CGAL::assign(e,result);
      f = get_visible_facet(e, Ray_3(p, s.target()));	
      if( f != Halffacet_handle())
	return f->incident_volume();
      SM_decorator SD(&*v); // now, the vertex has no incident facets
      CGAL_assertion( SD.number_of_sfaces() == 1);
      return SD.sfaces_begin()->volume();      
*/
    } else if( CGAL::assign( f, result)) {
      _CGAL_NEF_TRACEN("facet hit, obtaining volume...");
      if(f->plane().oriented_side(p) == ON_NEGATIVE_SIDE)
	f = f->twin();
      return Object_handle(f->incident_volume());
#ifdef CGAL_NEF3_FACET_WITH_BOX
    } else if( CGAL::assign(pf, *o)) {
      CGAL_assertion_msg(false, "should not be executed");
      Halffacet_handle f = pf.f;
      if(f->plane().oriented_side(p) == ON_NEGATIVE_SIDE)
	f = f->twin();
      return Object_handle(f->incident_volume()); 
#endif
#ifdef CGAL_NEF3_TRIANGULATE_FACETS
    } else if( CGAL::assign(t, result)) {
      f = t;
      _CGAL_NEF_TRACEN("facet hit, obtaining volume...");
      if(f->plane().oriented_side(p) == ON_NEGATIVE_SIDE)
	f = f->twin();
      return Object_handle(f->incident_volume());
#endif
    } else if( CGAL::assign(e, result)) {
      SM_decorator SD(&*e->source());
      if( SD.is_isolated(e))
        return Object_handle(e->incident_sface()->volume());	
      return get_visible_facet(e,Ray_3(s.source(),s.to_vector()))->incident_volume();
    }
    CGAL_assertion_msg(false, "wrong handle type");
    return Object_handle();
  }
  }