- Ifc: support IfcPolygonalBoundedHalfSpace entities (WIP).

code/IFCUtil.cpp

@@ -570,6 +570,7 @@ void ConvertTransformOperator(IfcMatrix4& out, const IfcCartesianTransformationO
 	out = locm * out * s;
 }
 
+
 } // ! IFC
 } // ! Assimp
 

code/IFCUtil.h

@@ -299,6 +299,7 @@ bool GenerateOpenings(std::vector<TempOpening>& openings,
 					  const IfcVector3& wall_extrusion_axis = IfcVector3(0,1,0));
 
 
+
 // IFCCurve.cpp
 
 // ------------------------------------------------------------------------------------------------
@@ -403,7 +404,6 @@ public:
 
 // IfcProfile.cpp
 bool ProcessCurve(const IfcCurve& curve,  TempMesh& meshout, ConversionData& conv);
-
 }
 }
 

code/IfcBoolean.cpp

@@ -182,6 +182,108 @@ void ProcessBooleanHalfSpaceDifference(const IfcHalfSpaceSolid* hs, TempMesh& re
 	IFCImporter::LogDebug("generating CSG geometry by plane clipping (IfcBooleanClippingResult)");
 }
 
+// ------------------------------------------------------------------------------------------------
+// Check if e0-e1 intersects a sub-segment of the given boundary line.
+// note: this method works on 3D vectors, but performs its intersection checks solely in xy.
+bool IntersectsBoundaryProfile( const IfcVector3& e0, const IfcVector3& e1, const std::vector<IfcVector3>& boundary,
+	std::vector<size_t>& intersected_boundary_segments,
+	std::vector<IfcVector3>& intersected_boundary_points,
+	bool half_open = false)
+{
+	ai_assert(intersected_boundary_segments.empty());
+	ai_assert(intersected_boundary_points.empty());
+
+	const IfcVector3& e = e1 - e0;
+
+	for (size_t i = 0, bcount = boundary.size(); i < bcount; ++i) {
+		// boundary segment i: b0-b1
+		const IfcVector3& b0 = boundary[i];
+		const IfcVector3& b1 = boundary[(i+1) % bcount];
+
+		const IfcVector3& b = b1 - b0;
+
+		// segment-segment intersection
+		// solve b0 + b*s = e0 + e*s for (s,t)
+		const IfcFloat det = (-b.x * e.y + e.x * b.y);
+		if(fabs(det) < 1e-6) {
+			// no solutions (parallel lines)
+			continue;
+		}
+
+		const IfcFloat x = b0.x - e0.x;
+		const IfcFloat y = b0.y - e0.y;
+
+		const IfcFloat s = (x*e.y - e.x*y)/det;
+		const IfcFloat t = (x*b.y - b.x*y)/det;
+
+#ifdef _DEBUG
+		const IfcVector3 check = b0 + b*s  - (e0 + e*t);
+		ai_assert((IfcVector2(check.x,check.y)).SquareLength() < 1e-5);
+#endif
+
+		// for a valid intersection, s-t should be in range [0,1]
+		if (s >= 0.0 && (s <= 1.0 || half_open) && t >= 0.0 && t <= 1.0) {
+	
+			const IfcVector3& p = b0 + b*s;
+
+			// only insert the point into the list if it is sufficiently
+			// far away from the previous intersection point. This way,
+			// we avoid duplicate detection if the intersection is
+			// directly on the vertex between two segments.
+			if (!intersected_boundary_points.empty() && intersected_boundary_segments.back()==(i==0?bcount-1:i-1) ) {
+				if((intersected_boundary_points.back() - p).SquareLength() < 1e-5) {
+					continue;
+				}
+			}
+			intersected_boundary_segments.push_back(i);
+			intersected_boundary_points.push_back(p);
+		}
+	}
+
+	return false;
+}
+
+
+// ------------------------------------------------------------------------------------------------
+bool PointInPoly(const IfcVector3& p, const std::vector<IfcVector3>& boundary)
+{
+	// even-odd algorithm: take a random vector that extends from p to infinite
+	// and counts how many times it intersects edges of the boundary.
+	// because checking for segment intersections is prone to numeric inaccuracies
+	// or double detections (i.e. when hitting multiple adjacent segments at their
+	// shared vertices) we do it thrice with different rays and vote on it.
+
+	std::vector<size_t> intersected_boundary_segments;
+	std::vector<IfcVector3> intersected_boundary_points;
+	size_t votes = 0;
+
+	IntersectsBoundaryProfile(p, p + IfcVector3(1.0,0,0), boundary, 
+		intersected_boundary_segments, 
+		intersected_boundary_points, true);
+
+	votes += intersected_boundary_segments.size() % 2;
+
+	intersected_boundary_segments.clear();
+	intersected_boundary_points.clear();
+
+	IntersectsBoundaryProfile(p, p + IfcVector3(0,1.0,0), boundary, 
+		intersected_boundary_segments, 
+		intersected_boundary_points, true);
+
+	votes += intersected_boundary_segments.size() % 2;
+
+	intersected_boundary_segments.clear();
+	intersected_boundary_points.clear();
+
+	IntersectsBoundaryProfile(p, p + IfcVector3(0,0,1.0), boundary, 
+		intersected_boundary_segments, 
+		intersected_boundary_points, true);
+
+	votes += intersected_boundary_segments.size() % 2;
+	return votes > 1;
+}
+
+
 // ------------------------------------------------------------------------------------------------
 void ProcessPolygonalBoundedBooleanHalfSpaceDifference(const IfcPolygonalBoundedHalfSpace* hs, TempMesh& result, 
 													   const TempMesh& first_operand, 
@@ -189,9 +291,249 @@ void ProcessPolygonalBoundedBooleanHalfSpaceDifference(const IfcPolygonalBounded
 {
 	ai_assert(hs != NULL);
 
-	return; // niy
+	const IfcPlane* const plane = hs->BaseSurface->ToPtr<IfcPlane>();
+	if(!plane) {
+		IFCImporter::LogError("expected IfcPlane as base surface for the IfcHalfSpaceSolid");
+		return;
+	}
+
+	// extract plane base position vector and normal vector
+	IfcVector3 p,n(0.f,0.f,1.f);
+	if (plane->Position->Axis) {
+		ConvertDirection(n,plane->Position->Axis.Get());
+	}
+	ConvertCartesianPoint(p,plane->Position->Location);
+
+	if(!IsTrue(hs->AgreementFlag)) {
+		n *= -1.f;
+	}
+
+	n.Normalize();
+
+	// obtain the polygonal bounding volume
+	boost::shared_ptr<TempMesh> profile = boost::shared_ptr<TempMesh>(new TempMesh());
+	if(!ProcessCurve(hs->PolygonalBoundary, *profile.get(), conv)) {
+		IFCImporter::LogError("expected valid polyline for boundary of boolean halfspace");
+		return;
+	}
+
+	IfcMatrix4 mat;
+	ConvertAxisPlacement(mat,hs->Position);
+	profile->Transform(mat);
+
+	// project the profile onto the plane (orthogonally along the plane normal)
+	IfcVector3 r;
+	bool have_r = false;
+	BOOST_FOREACH(IfcVector3& vec, profile->verts) {
+		vec = vec + ((p - vec) * n) * n;
+		ai_assert(fabs((vec-p) * n) < 1e-6);
+
+		if (!have_r && (vec-p).SquareLength() > 1e-8) {
+			r = vec-p;
+			have_r = true;
+		}
+	}
+
+	if (!have_r) {
+		IFCImporter::LogError("polyline for boundary of boolean halfspace is degenerate");
+		return;
+	}
+
+	// and map everything into a plane coordinate space so all intersection
+	// tests can be done in 2D space.
+	IfcMatrix4 proj;
+
+	r.Normalize();
+
+	IfcVector3 u = n ^ r;
+	u.Normalize();
+
+	proj.a1 = r.x;
+	proj.a2 = r.y;
+	proj.a3 = r.z;
+	proj.b1 = u.x;
+	proj.b2 = u.y;
+	proj.b3 = u.z;
+	proj.c1 = n.x;
+	proj.c2 = n.y;
+	proj.c3 = n.z;
+	
+	BOOST_FOREACH(IfcVector3& vec, profile->verts) {
+		vec *= proj;
+	}
+
+	const IfcMatrix4 proj_inv = IfcMatrix4(proj).Inverse();
+
+	// clip the current contents of `meshout` against the plane we obtained from the second operand
+	const std::vector<IfcVector3>& in = first_operand.verts;
+	std::vector<IfcVector3>& outvert = result.verts;
+
+	std::vector<unsigned int>::const_iterator begin = first_operand.vertcnt.begin(), 
+		end = first_operand.vertcnt.end(), iit;
+
+	outvert.reserve(in.size());
+	result.vertcnt.reserve(first_operand.vertcnt.size());
+
+	std::vector<size_t> intersected_boundary_segments;
+	std::vector<IfcVector3> intersected_boundary_points;
+
+	unsigned int vidx = 0;
+	for(iit = begin; iit != end; vidx += *iit++) {
+		if (!*iit) {
+			continue;
+		}
+
+		unsigned int newcount = 0;
+		bool was_outside_boundary = !PointInPoly(in[vidx], profile->verts);
+
+		size_t last_intersected_boundary_segment;
+		IfcVector3 last_intersected_boundary_point;
+
+		bool extra_point_flag = false;
+		IfcVector3 extra_point;
+
+		for(unsigned int i = 0; i < *iit; ++i) {
+			// current segment: [i,i+1 mod size] or [*extra_point,i] if extra_point_flag is set
+			const IfcVector3& e0 = extra_point_flag ? extra_point : in[vidx+i];
+			const IfcVector3& e1 = extra_point_flag ? in[vidx+i] : in[vidx+(i+1)%*iit];
 
+			// does the current segment intersect the polygonal boundary?
+			const IfcVector3& e0_plane = proj * e0;
+			const IfcVector3& e1_plane = proj * e1;
+
+			intersected_boundary_segments.clear();
+			intersected_boundary_points.clear();
+
+			const bool is_boundary_intersection = IntersectsBoundaryProfile(e0_plane, e1_plane, profile->verts, 
+				intersected_boundary_segments, 
+				intersected_boundary_points);
+
+			const bool is_outside_boundary = is_boundary_intersection ? !was_outside_boundary : was_outside_boundary;
+
+			// does the current segment intersect the plane?
+			// (no extra check if this is an extra point)
+			IfcVector3 isectpos;
+			const Intersect isect = extra_point_flag ? Intersect_No : IntersectSegmentPlane(p,n,e0,e1,isectpos);
+
+			// is it on the side of the plane that we keep?
+			const bool is_white_side =(e0-p).Normalize()*n > 0;
+
+			// e0 on good side of plane? (i.e. we should keep geometry on this side)
+			if (is_white_side) {
+				// but is there an intersection in e0-e1 and is e1 in the clipping
+				// boundary? In this case, generate a line that only goes to the
+				// intersection point.
+				if (isect == Intersect_Yes && PointInPoly(e1, profile->verts)) {
+					outvert.push_back(e0);
+					++newcount;
+
+					outvert.push_back(isectpos);
+					++newcount;
+
+					// this is, however, only a line that goes to the plane, but not
+					// necessarily to the point where the bounding volume on the
+					// black side of the plane is hit. So basically, we need another 
+					// check for [isectpos-e1], which should give an intersection
+					// point and also set the last_intersected_boundary_*'s.
+					extra_point_flag = true;
+					extra_point = isectpos;
+					continue;
+				}
+				else {
+					outvert.push_back(e0);
+					++newcount;
+				}
+			}
+			// e0 on bad side of plane (i.e. we should remove geometry on this side,
+			// but only if it is within the bounding volume).
+			else if (isect == Intersect_Yes) {
+				if (is_boundary_intersection) {}
+				// drop it and keep e1 instead
+				outvert.push_back(isectpos);
+				++newcount;
+			}
+			else {
+
+				// did we just pass the boundary line?
+				if (is_boundary_intersection) {
+
+					// and are now outside the clipping boundary?
+					if (is_outside_boundary) {
+						// in this case, get the point where the clipping boundary
+						// was entered first. Then, get the point where the clipping
+						// boundary volume was left! These two points with the plane
+						// normal form another plane that intersects the clipping
+						// volume. There are two ways to get from the first to the
+						// second point along the intersection curve, try to pick the
+						// one that lies within the current polygon.
+
+						// TODO this approach doesn't handle all cases
+
+						// ...
+
+						outvert.push_back(proj_inv * intersected_boundary_points.back());
+						++newcount;
+
+						outvert.push_back(e1);
+						++newcount;
+					}
+					else {
+						// we just entered the clipping boundary. Record the point
+						// and the segment where we entered and also generate this point.
+						last_intersected_boundary_segment = intersected_boundary_segments.front();
+						last_intersected_boundary_point = intersected_boundary_points.front();
+
+						outvert.push_back(e0);
+						++newcount;
+
+						outvert.push_back(proj_inv * last_intersected_boundary_point);
+						++newcount;
+					}
+				}				
+				// if not, we just keep the vertex
+				else if (is_outside_boundary) {
+					outvert.push_back(e0);
+					++newcount;
+				}
+			}
+
+			was_outside_boundary = is_outside_boundary;
+			extra_point_flag = false;
+		}
 	
+		if (!newcount) {
+			continue;
+		}
+
+		IfcVector3 vmin,vmax;
+		ArrayBounds(&*(outvert.end()-newcount),newcount,vmin,vmax);
+
+		// filter our IfcFloat points - those may happen if a point lies
+		// directly on the intersection line. However, due to IfcFloat
+		// precision a bitwise comparison is not feasible to detect
+		// this case.
+		const IfcFloat epsilon = (vmax-vmin).SquareLength() / 1e6f;
+		FuzzyVectorCompare fz(epsilon);
+
+		std::vector<IfcVector3>::iterator e = std::unique( outvert.end()-newcount, outvert.end(), fz );
+
+		if (e != outvert.end()) {
+			newcount -= static_cast<unsigned int>(std::distance(e,outvert.end()));
+			outvert.erase(e,outvert.end());
+		}
+		if (fz(*( outvert.end()-newcount),outvert.back())) {
+			outvert.pop_back();
+			--newcount;
+		}
+		if(newcount > 2) {
+			result.vertcnt.push_back(newcount);
+		}
+		else while(newcount-->0) {
+			result.verts.pop_back();
+		}
+
+	}
+	IFCImporter::LogDebug("generating CSG geometry by plane clipping with polygonal bounding (IfcBooleanClippingResult)");	
 }
 
 // ------------------------------------------------------------------------------------------------
@@ -283,6 +625,7 @@ void ProcessBoolean(const IfcBooleanResult& boolean, TempMesh& result, Conversio
 		}
 
 		if(hs) {
+
 			const IfcPolygonalBoundedHalfSpace* const hs_bounded = clip->SecondOperand->ResolveSelectPtr<IfcPolygonalBoundedHalfSpace>(conv.db);
 			if (hs_bounded) {
 				ProcessPolygonalBoundedBooleanHalfSpaceDifference(hs_bounded, result, first_operand, conv);