oomph::PolygonHelper Namespace Reference

Namespace with helper function for polygons. More...

Functions
bool	point_is_in_polygon (const Vector< double > &point, const Vector< Vector< double > > &polygon_vertex)

Detailed Description

Namespace with helper function for polygons.

///////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////

Function Documentation

point_is_in_polygon()

bool oomph::PolygonHelper::point_is_in_polygon	(	const Vector< double > &	point,
		const Vector< Vector< double > > &	polygon_vertex
	)

Check if point in inside polygon. Reference: http://paulbourke.net/geometry/insidepoly/

 {
 
  // Total number of vertices
  unsigned nvertex=polygon_vertex.size();
 
// faire: obviously polygons should ALWAYS be designed to be closed, so
// such testing should only need to be done once upon implementation of
// the polygon, not everytime we want to check number of intersections 
// hierher #ifdef PARANOID
  
  // Make sure the polygon closes exactly:
  double d0=polygon_vertex[nvertex-1][0]-polygon_vertex[0][0];
  double d1=polygon_vertex[nvertex-1][1]-polygon_vertex[0][1];
  if (sqrt(d0*d0+d1*d1)>0.0)
   {
    std::stringstream junk;
    junk << "First and last point of polygon don't coincide!\n"
         << "First point at: " 
         << polygon_vertex[0][0] << " " 
         << polygon_vertex[0][1] << "\n"
         << "Last point at: " 
         << polygon_vertex[nvertex-1][0] << " " 
         << polygon_vertex[nvertex-1][1] << "\n";
    throw OomphLibError(
     junk.str(),
     OOMPH_CURRENT_FUNCTION,
     OOMPH_EXCEPTION_LOCATION);
   }
 
// #endif
 
  // Counter for number of intersections
  unsigned intersect_counter=0;
  
  //Get first vertex
  Vector<double> p1=polygon_vertex[0];
  for (unsigned i=1;i<=nvertex;i++)
   {
    // Get second vertex by wrap-around
    Vector<double> p2 = polygon_vertex[i%nvertex];
    
    if (point[1] > std::min(p1[1],p2[1]))
     {
      if (point[1] <= std::max(p1[1],p2[1]))
       {
        if (point[0] <= std::max(p1[0],p2[0]))
         {
          if (p1[1] != p2[1])
           {
            double xintersect =
             (point[1]-p1[1])*(p2[0]-p1[0])/
             (p2[1]-p1[1])+p1[0];
            if ( (p1[0] == p2[0]) ||
                 (point[0] <= xintersect) )
             {
              intersect_counter++;
             }
           }
         }
       }
     }
    p1 = p2;
   }
  
  // Even number of intersections: outside
  if (intersect_counter%2==0)
   {
    return false;
   }
  return true;
 }

References i, max, min, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, p1, and sqrt().