Find out if point is inside solid

[TIG]

First check if the point is within the solid.bounds...

If not, then it is not inside, although if it's inside then of course it might be on its 'surface' and not actually inside it !

If it is a true manifold solid then all of its face.normal's should point 'outwards'.

Use this test on the ' point':
Choose an arbitrary ' vector' like [1,0,0]
hit = model.raytest(point, vector)
If it returns no hit at all then it definitely is outside of everything - at least using that 'vector'!
If it returns a hit, then using that ' hit' array you can test if
hit.include?(solid)
If no hit then perhaps 'reverse' the vector and retest...
Assuming we get a hit that includes the 'solid'...
Check if the hit includes a 'face'.
If it's an edge, perhaps adjust the vector slightly to see if you can actually hit a face !
Assume you have hit a face, now compare the vector and the face.normal using
angle = vector.angle_between(face.normal)
If angle < 90.degrees then it is 'inside'.
If its >= then its outside.

[jiminy-billy-bob]

Maybe this could help: http://clb.demon.fi/MathGeoLib/nightly/docs/Polyhedron.cpp_code.html#532

@unknownuser said:

General strategy: pick a ray from the query point to a random direction, and count the number of times the ray intersects a face. If it intersects an odd number of times, the given point must have been inside the polyhedron.

[Garry K]

Randolf Franklin solved this problem back in the 1970's

Here is his original 2D code for a single polygon

int pnpoly(int nvert, float *vertx, float *verty, float testx, float testy)
{
  int i, j, c = 0;
  for (i = 0, j = nvert-1; i < nvert; j = i++) {
    if ( ((verty[i]>testy) != (verty[j]>testy)) &&
	 (testx < (vertx[j]-vertx[i]) * (testy-verty[i]) / (verty[j]-verty[i]) + vertx[i]) )
       c = !c;
  }
  return c;
}

I ran a modified version of this code on very old mobile computers 80186 processor running at 18 megahertz. We used it to determine which State or Province your vehicle is in.

3D isn't much more difficult but you do have to look at each face and see if your ray passes through it or not.

[sdmitch]

@icatic said:

Hello,

Is there a way using the ruby API to find out if a point is inside the actual volume of a solid instance/group?

Thanks,
Valdemar

Here is my solution. If it can't find an opening is 1000 tries, the point is most likely inside the solid.

def is_point_inside(pt,cg)
	return false unless cg.manifold?
	return false unless cg.bounds.contains?(pt)
	1000.times{
		vec=Geom;;Vector3d.new(rand*(rand<=>0.5),rand*(rand<=>0.5),rand*(rand<=>0.5))
		h=Sketchup.active_model.raytest([pt,vec]);
		return false unless (h && h[1].include?(cg))
	}
	return true
end

[fredo6]

Assuming you have a true solid,

• You have to test first if the point is on a face of the solid (you decide if it is 'in' or 'out')
• otherwise, you take one arbitrary oriented direction (say X direction) and check the number of times it crosses a face of the solid. It must be an odd number if inside, or even if outside. The trick is that you must keep track of the intersection points to make sure you do not reuse the same point twice or more times (for instance when the half-line crosses an edge or a vertex).

Fredo

[icatic]

Thanks for all the replies!

I ended up using the method that both Ruts and fredo6 explained. There's a special circumstance when the arbitrary line runs along a face parallell to the line that I had to consider in the code (the skip_intersection part). Here's the code I ended up with:

def solid_contains_point?(solid, point, include_shell = true)
  return false unless solid.manifold?
  return false unless solid.bounds.contains? point
  point = point.transform(solid.transformation.inverse)
  line_vector = Geom;;Vector3d.new(1, 0, 0)
  line = [point, line_vector]
  previous_intersection_point = previous_classify_result = nil
  intersection_points = Array.new

  solid.definition.entities.find_all { |e| e.instance_of? Sketchup;;Face }.each { |face|

    intersection_point = Geom.intersect_line_plane(line, face.plane)
    if intersection_point and intersection_point != previous_intersection_point
      classify_result = face.classify_point(intersection_point)
      point_on_face = [Sketchup;;Face;;PointInside, Sketchup;;Face;;PointOnVertex, Sketchup;;Face;;PointOnEdge].include? classify_result
      angle = line_vector.dot(point.vector_to(intersection_point))
      if point_on_face and ((angle >= 0 and include_shell) or (angle > 0 and not include_shell))
        intersection_points.push [intersection_point, classify_result]
        previous_intersection_point = intersection_point
      end
    end

  }

  previous_skip_intersection = false
  intersections = 0
  intersection_points.sort! { |a,b| point.vector_to(a[0]).length <=> point.vector_to(b[0]).length }.each { |array|

    intersection_point = array[0]
    classify_result = array[1]

    # skip intersection point if there's another face that contains both the current
    # intersection point and the previous intersection point
    skip_intersection = false
    if previous_classify_result and previous_classify_result != Sketchup;;Face;;PointInside
        solid.definition.entities.find_all { |e| e.instance_of? Sketchup;;Face }.each { |e|
          if classify_result == Sketchup;;Face;;PointOnVertex
            if e.vertices.find_all { |vertex| vertex.position == intersection_point or vertex.position == previous_intersection_point }.size == 2
              skip_intersection = true
              break
            end
          elsif classify_result == Sketchup;;Face;;PointOnEdge
            if e.edges.find_all { |edge| intersection_point.on_line?(edge.line) or previous_intersection_point.on_line?(edge.line) }.size == 2
              skip_intersection = true
              break
            end
          end
        }
    end

    # count all intersection points connected by face parallell to and on line as one intersection
    intersections += 1 if previous_skip_intersection and !skip_intersection
    
    previous_skip_intersection = skip_intersection
    next if skip_intersection

    intersections += 1
    previous_intersection_point = intersection_point
    previous_classify_result = classify_result

  }

  intersections.odd?
end

[fredo6]

That should make it.

Just 2 remarks

• I would cache the list of valid classify_status outside the loop
• I would construct the list of intersection points without bothering whether this is a vertex or on edge and then sort the list by the distance to the point or X, Y, Z coordinates and eliminate the redundant points (which are necessarily in sequence)

Fredo

[Paciente8159]

Is this OK?
The same principle. A raycast from the testing point. I'm only counting how many faces it sees in one direction and only if the point hits the face on the inside (from the testing point)
To also include the outter shell of the solid the point classification must be changed to validade Vertex and Edges hits.

Thanks

Here is the code

def point_insidesolid?(point, solid)
		return false unless solid.manifold?

		#vector that goes through the origin
		vec = Geom;;Vector3d.new(1,0,0).normalize
		counter = 0
		solid.entities.grep(Sketchup;;Face).each do |face|
			pt = Geom.intersect_line_plane([point,vec], [face.vertices[0].position, face.normal])
			if(!pt.nil?)
				if(face.classify_point(pt)==Sketchup;;Face;;PointInside)
					if(pt.x>point.x)
						counter +=1
					end
				end
			end
		end
		
		return false unless counter%2!=0
		return true
	end

[sdmitch]

@paciente8159 said:

Is this OK?
The same principle. A raycast from the testing point. I'm only counting how many faces it sees in one direction and only if the point hits the face on the inside (from the testing point)
To also include the outter shell of the solid the point classification must be changed to validade Vertex and Edges hits.

Thanks

Here is the code

def point_insidesolid?(point, solid)
> 		return false unless solid.manifold?
> 
> 		#vector that goes through the origin
> 		vec = Geom;;Vector3d.new(1,0,0).normalize
> 		counter = 0
> 		solid.entities.grep(Sketchup;;Face).each do |face|
> 			pt = Geom.intersect_line_plane([point,vec], [face.vertices[0].position, face.normal])
> 			if(!pt.nil?)
> 				if(face.classify_point(pt)==Sketchup;;Face;;PointInside)
> 					if(pt.x>point.x)
> 						counter +=1
> 					end
> 				end
> 			end
> 		end
> 		
> 		return false unless counter%2!=0
> 		return true
> 	end

It should if "point" is truly inside the solid. If it is a model point then it will need to be transformed by the inverse of the solid's transformation.

[Paciente8159]

@sdmitch said:

It should if "point" is truly inside the solid. If it is a model point then it will need to be transformed by the inverse of the solid's transformation.

Thanks
Yes the point must be in the solid(group) coordinate space.

I've changed the code a bit to allow shell (outer face) testing
Here it is:

Edited:
After some testing still something was not OK.
I realized that when the ray crosses an edge it hits 2 different planes and through a vertice it can hit an unknown number of faces. The solution I came is that the ray crosses the solid in only one point. So instead of keeping track of how many faces it hits, I register all points the ray hits the solid and then perform a uniq-like operation to the array to count the number of points/hits.
I've updated the code.

For the filtering uniq! doesn't work (guessing because it uses eql? and this is not implemented in Point3d)
Also tried

pts.uniq! {|pt| pt.to_a}

but it didn't worked, but since the raycast is parallel to X axis comparing x coordinate does the trick

def point_insidesolid?(point, solid, exclude_surface=false)		
		vect = Geom;;Vector3d.new(1,0,0)
		pts = []
		solid.entities.grep(Sketchup;;Face).each do |face|
			if(FaceTools.is_pointof_face?(point,face))
				return true unless exclude_surface
				return false
			else
				pt = Geom.intersect_line_plane([point,vect], [face.vertices[0].position, face.normal])
				if(!pt.nil?)
					if(FaceTools.is_pointof_face?(pt,face) && pt.x>point.x)
						pts.push pt
					end
				end
			end
		end
		
		ptsuniq = pts.uniq {|pt| pt.x}

		return ptsuniq.count%2!=0
	end

	def point_insideface?(point,face)
		case(face.classify_point(point))
			when Sketchup;;Face;;PointInside, Sketchup;;Face;;PointOnVertex, Sketchup;;Face;;PointOnEdge
				return true
		end
		
		return false
	end