Thread: [brlcad-commits] SF.net SVN: brlcad:[55433] brlcad/trunk/src/libbrep/intersect.cpp
Open Source Solid Modeling CAD
Brought to you by:
brlcad
From: <pho...@us...> - 2013-05-15 17:55:10
|
Revision: 55433 http://sourceforge.net/p/brlcad/code/55433 Author: phoenixyjll Date: 2013-05-15 17:55:05 +0000 (Wed, 15 May 2013) Log Message: ----------- Use max_dis_2dA and max_dis_2dB for merging polylines. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-05-15 17:38:44 UTC (rev 55432) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-05-15 17:55:05 UTC (rev 55433) @@ -192,9 +192,9 @@ struct PointPair { int indexA, indexB; - double distance; + double distance3d, distanceA, distanceB; bool operator < (const PointPair &_pp) const { - return distance < _pp.distance; + return distance3d < _pp.distance3d; } }; @@ -504,8 +504,8 @@ max_dis = pow(surfA->BoundingBox().Volume()*surfB->BoundingBox().Volume(), 1.0/6.0) * 0.2; } - double max_dis_2dA = ON_2dVector(surfA->Domain(0).Length(), surfA->Domain(1).Length()).Length() * 0.01; - double max_dis_2dB = ON_2dVector(surfB->Domain(0).Length(), surfB->Domain(1).Length()).Length() * 0.01; + double max_dis_2dA = ON_2dVector(surfA->Domain(0).Length(), surfA->Domain(1).Length()).Length() * 0.1; + double max_dis_2dB = ON_2dVector(surfB->Domain(0).Length(), surfB->Domain(1).Length()).Length() * 0.1; bu_log("max_dis: %lf\n", max_dis); bu_log("max_dis_2dA: %lf\n", max_dis_2dA); bu_log("max_dis_2dB: %lf\n", max_dis_2dB); @@ -515,8 +515,10 @@ for (int i = 0; i < curvept.Count(); i++) { for (int j = i + 1; j < curvept.Count(); j++) { PointPair ppair; - ppair.distance = curvept[i].DistanceTo(curvept[j]); - if (ppair.distance < max_dis) { + ppair.distance3d = curvept[i].DistanceTo(curvept[j]); + ppair.distanceA = curveuv[i].DistanceTo(curveuv[j]); + ppair.distanceB = curvest[i].DistanceTo(curvest[j]); + if (ppair.distanceA < max_dis_2dA && ppair.distanceB < max_dis_2dB) { ppair.indexA = i; ppair.indexB = j; ptpairs.push_back(ppair); @@ -682,3 +684,12 @@ } #endif + +// Local Variables: +// tab-width: 8 +// mode: C++ +// c-basic-offset: 4 +// indent-tabs-mode: t +// c-file-style: "stroustrup" +// End: +// ex: shiftwidth=4 tabstop=8 This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-03 13:48:48
|
Revision: 55637 http://sourceforge.net/p/brlcad/code/55637 Author: phoenixyjll Date: 2013-06-03 13:48:45 +0000 (Mon, 03 Jun 2013) Log Message: ----------- remove the unused variable names to quell compiler warnings. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-03 13:43:43 UTC (rev 55636) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-03 13:48:45 UTC (rev 55637) @@ -64,11 +64,11 @@ * Point-curve intersections (PCI) */ bool -ON_Intersect(const ON_3dPoint& pointA, - const ON_Curve& curveB, - ON_ClassArray<ON_PX_EVENT>& x, - double tolerance, - const ON_Interval* curveB_domain) +ON_Intersect(const ON_3dPoint&, + const ON_Curve&, + ON_ClassArray<ON_PX_EVENT>&, + double, + const ON_Interval*) { // Implement later. return false; @@ -78,12 +78,12 @@ * Point-surface intersections (PSI) */ bool -ON_Intersect(const ON_3dPoint& pointA, - const ON_Surface& surfaceB, - ON_ClassArray<ON_PX_EVENT>& x, - double tolerance, - const ON_Interval* surfaceB_udomain, - const ON_Interval* surfaceB_vdomain) +ON_Intersect(const ON_3dPoint&, + const ON_Surface&, + ON_ClassArray<ON_PX_EVENT>&, + double, + const ON_Interval*, + const ON_Interval*) { // Implement later. return false; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <br...@us...> - 2013-06-06 15:59:21
|
Revision: 55682 http://sourceforge.net/p/brlcad/code/55682 Author: brlcad Date: 2013-06-06 15:59:17 +0000 (Thu, 06 Jun 2013) Log Message: ----------- bio.h is a system header wrapper but should still probably come after the <> headers Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-05 22:32:42 UTC (rev 55681) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-06 15:59:17 UTC (rev 55682) @@ -25,16 +25,17 @@ #include "common.h" -#include "bio.h" #include <assert.h> #include <vector> #include <algorithm> +#include "bio.h" #include "vmath.h" #include "bu.h" #include "brep.h" + /** * Point-point intersections (PPI) */ This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-07 16:59:34
|
Revision: 55692 http://sourceforge.net/p/brlcad/code/55692 Author: phoenixyjll Date: 2013-06-07 16:59:29 +0000 (Fri, 07 Jun 2013) Log Message: ----------- Try to add PS support using get_closest_point(). Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-07 15:29:06 UTC (rev 55691) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-07 16:59:29 UTC (rev 55692) @@ -38,6 +38,14 @@ /** * Point-point intersections (PPI) + * + * approach: + * + * - Calculate the distance of the two points. + * + * - If the distance is within the tolerance, the two points + * intersect. Thd mid-point of them is the intersection point, + * and half of the distance is the uncertainty radius. */ bool ON_Intersect(const ON_3dPoint& pointA, @@ -79,14 +87,49 @@ * Point-surface intersections (PSI) */ bool -ON_Intersect(const ON_3dPoint&, - const ON_Surface&, - ON_ClassArray<ON_PX_EVENT>&, - double, - const ON_Interval*, - const ON_Interval*) +ON_Intersect(const ON_3dPoint& pointA, + const ON_Surface& surfaceB, + ON_ClassArray<ON_PX_EVENT>& x, + double tolerance, + const ON_Interval* surfaceB_udomain, + const ON_Interval* surfaceB_vdomain) { - // Implement later. + if (tolerance <= 0.0) + tolerance = 0.01; + + ON_Interval u_domain, v_domain; + if (surfaceB_udomain == 0) + u_domain = surfaceB.Domain(0); + else + u_domain = *surfaceB_udomain; + if (surfaceB_vdomain == 0) + v_domain = surfaceB.Domain(1); + else + v_domain = *surfaceB_vdomain; + + // We need ON_BrepFace for get_closest_point(). + // TODO: Use routines like SubSurface in SSI to reduce computation. + ON_Brep *brep = ON_Brep::New(); + brep->AddSurface(surfaceB.Duplicate()); + brep->NewFace(0); + ON_2dPoint closest_point_uv; + brlcad::get_closest_point(closest_point_uv, brep->Face(0), pointA); + delete brep; + + ON_3dPoint closest_point_3d = surfaceB.PointAt(closest_point_uv.x, closest_point_uv.y); + bu_log("%lf %lf %lf\n", pointA.x, pointA.y, pointA.z); + bu_log("%lf %lf %lf\n", closest_point_3d.x, closest_point_3d.y, closest_point_3d.z); + if (pointA.DistanceTo(closest_point_3d) <= tolerance) { + ON_PX_EVENT Event; + Event.m_type = ON_PX_EVENT::psx_point; + Event.m_A = pointA; + Event.m_b = closest_point_uv; + Event.m_B = closest_point_3d; + Event.m_Mid = (pointA + Event.m_B) * 0.5; + Event.m_radius = pointA.DistanceTo(Event.m_B) * 0.5; + x.Append(Event); + return true; + } return false; } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-12 07:12:54
|
Revision: 55715 http://sourceforge.net/p/brlcad/code/55715 Author: phoenixyjll Date: 2013-06-12 07:12:47 +0000 (Wed, 12 Jun 2013) Log Message: ----------- Calculate point-curve intersection using sub-division. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-11 21:55:33 UTC (rev 55714) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-12 07:12:47 UTC (rev 55715) @@ -36,6 +36,168 @@ #include "brep.h" +/* Sub-division support for a curve. + * It's similar to generating the bounding box tree, when the Split() + * method is called, the curve is splitted into two parts, whose bounding + * boxes become the children of the original curve's bbox. + */ +struct Subcurve { +private: + ON_BoundingBox m_node; +public: + ON_Curve *m_curve; + ON_Interval m_t; + Subcurve *m_children[2]; + ON_BOOL32 m_islinear; + + Subcurve() : m_curve(NULL), m_islinear(false) + { + m_children[0] = m_children[1] = NULL; + } + Subcurve(const Subcurve &_scurve) + { + m_islinear = _scurve.m_islinear; + m_curve = _scurve.m_curve->Duplicate(); + m_t = _scurve.m_t; + m_children[0] = m_children[1] = NULL; + SetBBox(_scurve.m_node); + } + ~Subcurve() + { + for (int i = 0; i < 2; i++) { + if (m_children[i]) + delete m_children[i]; + } + delete m_curve; + } + int Split() + { + for (int i = 0; i < 2; i++) + m_children[i] = new Subcurve(); + + if (m_curve->Split(m_t.Mid(), m_children[0]->m_curve, m_children[1]->m_curve) == false) + return -1; + + m_children[0]->m_t = ON_Interval(m_t.Min(), m_t.Mid()); + m_children[0]->SetBBox(m_children[0]->m_curve->BoundingBox()); + m_children[0]->m_islinear = m_children[0]->m_curve->IsLinear(); + m_children[1]->m_t = ON_Interval(m_t.Mid(), m_t.Max()); + m_children[1]->SetBBox(m_children[1]->m_curve->BoundingBox()); + m_children[1]->m_islinear = m_children[1]->m_curve->IsLinear(); + + return 0; + } + void GetBBox(ON_3dPoint &min, ON_3dPoint &max) + { + min = m_node.m_min; + max = m_node.m_max; + } + void SetBBox(const ON_BoundingBox &bbox) + { + m_node = bbox; + } + bool IsPointIn(const ON_3dPoint &pt) + { + return m_node.IsPointIn(pt); + } +}; + + +/* Sub-division support for a surface. + * It's similar to generating the bounding box tree, when the Split() + * method is called, the surface is splitted into two parts, whose bounding + * boxes become the children of the original surface's bbox. + */ +struct Subsurface { +private: + ON_BoundingBox m_node; +public: + ON_Surface *m_surf; + ON_Interval m_u, m_v; + Subsurface *m_children[4]; + + Subsurface() : m_surf(NULL) + { + m_children[0] = m_children[1] = m_children[2] = m_children[3] = NULL; + } + Subsurface(const Subsurface &_ssurf) + { + m_surf = _ssurf.m_surf->Duplicate(); + m_u = _ssurf.m_u; + m_v = _ssurf.m_v; + m_children[0] = m_children[1] = m_children[2] = m_children[3] = NULL; + SetBBox(_ssurf.m_node); + } + ~Subsurface() + { + for (int i = 0; i < 4; i++) { + if (m_children[i]) + delete m_children[i]; + } + delete m_surf; + } + int Split() + { + for (int i = 0; i < 4; i++) + m_children[i] = new Subsurface(); + ON_Surface *temp_surf1 = NULL, *temp_surf2 = NULL; + ON_BOOL32 ret = true; + ret = m_surf->Split(0, m_u.Mid(), temp_surf1, temp_surf2); + if (!ret) { + delete temp_surf1; + delete temp_surf2; + return -1; + } + + ret = temp_surf1->Split(1, m_v.Mid(), m_children[0]->m_surf, m_children[1]->m_surf); + delete temp_surf1; + if (!ret) { + delete temp_surf2; + return -1; + } + m_children[0]->m_u = ON_Interval(m_u.Min(), m_u.Mid()); + m_children[0]->m_v = ON_Interval(m_v.Min(), m_v.Mid()); + m_children[0]->SetBBox(m_children[0]->m_surf->BoundingBox()); + m_children[1]->m_u = ON_Interval(m_u.Min(), m_u.Mid()); + m_children[1]->m_v = ON_Interval(m_v.Mid(), m_v.Max()); + m_children[1]->SetBBox(m_children[1]->m_surf->BoundingBox()); + + ret = temp_surf2->Split(1, m_v.Mid(), m_children[2]->m_surf, m_children[3]->m_surf); + delete temp_surf2; + if (!ret) + return -1; + m_children[2]->m_u = ON_Interval(m_u.Mid(), m_u.Max()); + m_children[2]->m_v = ON_Interval(m_v.Min(), m_v.Mid()); + m_children[2]->SetBBox(m_children[2]->m_surf->BoundingBox()); + m_children[3]->m_u = ON_Interval(m_u.Mid(), m_u.Max()); + m_children[3]->m_v = ON_Interval(m_v.Mid(), m_v.Max()); + m_children[3]->SetBBox(m_children[3]->m_surf->BoundingBox()); + + return 0; + } + void GetBBox(ON_3dPoint &min, ON_3dPoint &max) + { + min = m_node.m_min; + max = m_node.m_max; + } + void SetBBox(const ON_BoundingBox &bbox) + { + m_node = bbox; + /* Make sure that each dimension of the bounding box is greater than + * ON_ZERO_TOLERANCE. + * It does the same work as building the surface tree when ray tracing + */ + for (int i = 0; i < 3; i++) { + double d = m_node.m_max[i] - m_node.m_min[i]; + if (ON_NearZero(d, ON_ZERO_TOLERANCE)) { + m_node.m_min[i] -= 0.001; + m_node.m_max[i] += 0.001; + } + } + } +}; + + /** * Point-point intersections (PPI) * @@ -69,20 +231,100 @@ return false; } + +#define MAX_PCI_DEPTH 8 /** * Point-curve intersections (PCI) */ bool -ON_Intersect(const ON_3dPoint&, - const ON_Curve&, - ON_ClassArray<ON_PX_EVENT>&, - double, - const ON_Interval*) +ON_Intersect(const ON_3dPoint& pointA, + const ON_Curve& curveB, + ON_ClassArray<ON_PX_EVENT>& x, + double tolerance, + const ON_Interval* curveB_domain) { - // Implement later. + bu_log("PCI called.\n"); + Subcurve root; + if (curveB_domain == NULL) { + root.m_curve = curveB.Duplicate(); + root.m_t = curveB.Domain(); + } + else { + // Use ON_Curve::Split() to get the sub-curve inside curveB_domain + ON_Curve *temp_curve1 = NULL, *temp_curve2 = NULL, *temp_curve3 = NULL; + if (!curveB.Split(curveB_domain->Min(), temp_curve1, temp_curve2)) + return false; + delete temp_curve1; + temp_curve1 = NULL; + if (!temp_curve2->Split(curveB_domain->Max(), temp_curve1, temp_curve3)) + return false; + delete temp_curve1; + delete temp_curve2; + root.m_curve = temp_curve3; + root.m_t = *curveB_domain; + } + + root.SetBBox(root.m_curve->BoundingBox()); + root.m_islinear = root.m_curve->IsLinear(); + + if (!root.IsPointIn(pointA)) + return false; + + std::vector<Subcurve*> candidates, next_candidates; + candidates.push_back(&root); + + // Find the sub-curves that are possibly intersected with the point. + for (int i = 0; i < MAX_PCI_DEPTH; i++) { + next_candidates.clear(); + for (unsigned int j = 0; j < candidates.size(); j++) { + if (candidates[j]->m_islinear) { + next_candidates.push_back(candidates[j]); + } else { + if (candidates[j]->Split() == 0) { + if (candidates[j]->m_children[0]->IsPointIn(pointA)) + next_candidates.push_back(candidates[j]->m_children[0]); + if (candidates[j]->m_children[1]->IsPointIn(pointA)) + next_candidates.push_back(candidates[j]->m_children[1]); + } else + next_candidates.push_back(candidates[j]); + } + } + candidates = next_candidates; + } + + for (unsigned int i = 0; i < candidates.size(); i++) { + // Use linear approximation to get the intersection point + ON_Line line(candidates[i]->m_curve->PointAtStart(), candidates[i]->m_curve->PointAtEnd()); + double t, dis; + line.ClosestPointTo(pointA, &t); + // Get the closest point on the line *segment* to the point + double closest_point_t; + if (t < 0) + closest_point_t = 0; + else if (t > 1) + closest_point_t = 1; + else + closest_point_t = t; + ON_3dPoint closest_point = line.PointAt(closest_point_t); + dis = pointA.DistanceTo(closest_point); + + if (dis <= tolerance) { + ON_PX_EVENT *Event = new ON_PX_EVENT; + Event->m_type = ON_PX_EVENT::pcx_point; + Event->m_A = pointA; + Event->m_B = closest_point; + Event->m_b[0] = candidates[i]->m_t.Min()+candidates[i]->m_t.Length()*closest_point_t; + Event->m_Mid = (pointA + Event->m_B) * 0.5; + Event->m_radius = dis * 0.5; + x.Append(*Event); + return true; + } + } + // All candidates have no intersection return false; } + /** * Point-surface intersections (PSI) */ @@ -296,97 +538,7 @@ }; -struct Subsurface { -private: - ON_BoundingBox m_node; -public: - ON_Surface *m_surf; - ON_Interval m_u, m_v; - Subsurface *m_children[4]; - - Subsurface() : m_surf(NULL) - { - m_children[0] = m_children[1] = m_children[2] = m_children[3] = NULL; - } - Subsurface(const Subsurface &_ssurf) - { - m_surf = _ssurf.m_surf->Duplicate(); - m_u = _ssurf.m_u; - m_v = _ssurf.m_v; - m_children[0] = m_children[1] = m_children[2] = m_children[3] = NULL; - SetBBox(_ssurf.m_node); - } - ~Subsurface() - { - for (int i = 0; i < 4; i++) { - if (m_children[i]) - delete m_children[i]; - } - delete m_surf; - } - int Split() - { - for (int i = 0; i < 4; i++) - m_children[i] = new Subsurface(); - ON_Surface *temp_surf1 = NULL, *temp_surf2 = NULL; - ON_BOOL32 ret = true; - ret = m_surf->Split(0, m_u.Mid(), temp_surf1, temp_surf2); - if (!ret) { - delete temp_surf1; - delete temp_surf2; - return -1; - } - - ret = temp_surf1->Split(1, m_v.Mid(), m_children[0]->m_surf, m_children[1]->m_surf); - delete temp_surf1; - if (!ret) { - delete temp_surf2; - return -1; - } - m_children[0]->m_u = ON_Interval(m_u.Min(), m_u.Mid()); - m_children[0]->m_v = ON_Interval(m_v.Min(), m_v.Mid()); - m_children[0]->SetBBox(m_children[0]->m_surf->BoundingBox()); - m_children[1]->m_u = ON_Interval(m_u.Min(), m_u.Mid()); - m_children[1]->m_v = ON_Interval(m_v.Mid(), m_v.Max()); - m_children[1]->SetBBox(m_children[1]->m_surf->BoundingBox()); - - ret = temp_surf2->Split(1, m_v.Mid(), m_children[2]->m_surf, m_children[3]->m_surf); - delete temp_surf2; - if (!ret) - return -1; - m_children[2]->m_u = ON_Interval(m_u.Mid(), m_u.Max()); - m_children[2]->m_v = ON_Interval(m_v.Min(), m_v.Mid()); - m_children[2]->SetBBox(m_children[2]->m_surf->BoundingBox()); - m_children[3]->m_u = ON_Interval(m_u.Mid(), m_u.Max()); - m_children[3]->m_v = ON_Interval(m_v.Mid(), m_v.Max()); - m_children[3]->SetBBox(m_children[3]->m_surf->BoundingBox()); - - return 0; - } - void GetBBox(ON_3dPoint &min, ON_3dPoint &max) - { - min = m_node.m_min; - max = m_node.m_max; - } - void SetBBox(const ON_BoundingBox &bbox) - { - m_node = bbox; - /* Make sure that each dimension of the bounding box is greater than - * ON_ZERO_TOLERANCE. - * It does the same work as building the surface tree when ray tracing - */ - for (int i = 0; i < 3; i++) { - double d = m_node.m_max[i] - m_node.m_min[i]; - if (ON_NearZero(d, ON_ZERO_TOLERANCE)) { - m_node.m_min[i] -= 0.001; - m_node.m_max[i] += 0.001; - } - } - } -}; - - -#define INTERSECT_MAX_DEPTH 8 +#define MAX_SSI_DEPTH 8 int ON_Intersect(const ON_Surface* surfA, const ON_Surface* surfB, @@ -429,7 +581,7 @@ * method of splitting a NURBS surface. * So finally only a small subset of the surface tree is created. */ - for (int h = 0; h <= INTERSECT_MAX_DEPTH; h++) { + for (int h = 0; h <= MAX_SSI_DEPTH; h++) { if (nodepairs.empty()) break; NodePairs tmp_pairs; @@ -444,7 +596,7 @@ if (ret1) { if (ret2) { /* both splits failed */ tmp_pairs.push_back(*i); - h = INTERSECT_MAX_DEPTH; + h = MAX_SSI_DEPTH; } else { /* the first failed */ for (int j = 0; j < 4; j++) tmp_pairs.push_back(std::make_pair((*i).first, (*i).second->m_children[j])); @@ -469,7 +621,7 @@ (*i).first->GetBBox(box_a.m_min, box_a.m_max); (*i).second->GetBBox(box_b.m_min, box_b.m_max); if (box_intersect.Intersection(box_a, box_b)) { - if (h == INTERSECT_MAX_DEPTH) { + if (h == MAX_SSI_DEPTH) { // We have arrived at the bottom of the trees. // Get an estimate of the intersection point lying on the intersection curve This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-12 07:16:17
|
Revision: 55716 http://sourceforge.net/p/brlcad/code/55716 Author: phoenixyjll Date: 2013-06-12 07:16:14 +0000 (Wed, 12 Jun 2013) Log Message: ----------- ws. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-12 07:12:47 UTC (rev 55715) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-12 07:16:14 UTC (rev 55716) @@ -307,7 +307,7 @@ closest_point_t = t; ON_3dPoint closest_point = line.PointAt(closest_point_t); dis = pointA.DistanceTo(closest_point); - + if (dis <= tolerance) { ON_PX_EVENT *Event = new ON_PX_EVENT; Event->m_type = ON_PX_EVENT::pcx_point; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-17 14:38:56
|
Revision: 55793 http://sourceforge.net/p/brlcad/code/55793 Author: phoenixyjll Date: 2013-06-17 14:38:53 +0000 (Mon, 17 Jun 2013) Log Message: ----------- The tolerance should be considered in the IsPointIn() test. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-15 17:28:43 UTC (rev 55792) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-17 14:38:53 UTC (rev 55793) @@ -96,9 +96,11 @@ { m_node = bbox; } - bool IsPointIn(const ON_3dPoint &pt) + bool IsPointIn(const ON_3dPoint &pt, double tolerance = 0.0) { - return m_node.IsPointIn(pt); + ON_3dVector vtol(tolerance,tolerance,tolerance); + ON_BoundingBox new_bbox(m_node.m_min-vtol, m_node.m_max+vtol); + return new_bbox.IsPointIn(pt); } }; @@ -243,7 +245,9 @@ double tolerance, const ON_Interval* curveB_domain) { - bu_log("PCI called.\n"); + if (tolerance <= 0.0) + tolerance = 0.01; + Subcurve root; if (curveB_domain == NULL) { root.m_curve = curveB.Duplicate(); @@ -267,7 +271,7 @@ root.SetBBox(root.m_curve->BoundingBox()); root.m_islinear = root.m_curve->IsLinear(); - if (!root.IsPointIn(pointA)) + if (!root.IsPointIn(pointA, tolerance)) return false; std::vector<Subcurve*> candidates, next_candidates; @@ -281,9 +285,9 @@ next_candidates.push_back(candidates[j]); } else { if (candidates[j]->Split() == 0) { - if (candidates[j]->m_children[0]->IsPointIn(pointA)) + if (candidates[j]->m_children[0]->IsPointIn(pointA, tolerance)) next_candidates.push_back(candidates[j]->m_children[0]); - if (candidates[j]->m_children[1]->IsPointIn(pointA)) + if (candidates[j]->m_children[1]->IsPointIn(pointA, tolerance)) next_candidates.push_back(candidates[j]->m_children[1]); } else next_candidates.push_back(candidates[j]); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-18 08:40:46
|
Revision: 55795 http://sourceforge.net/p/brlcad/code/55795 Author: phoenixyjll Date: 2013-06-18 08:40:42 +0000 (Tue, 18 Jun 2013) Log Message: ----------- Use macros to represent default tolerance and change it to 0.001 (the same as the default tolerance of curve/curve, curve/surface, surface/surface defined by openNURBS). Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-17 14:43:25 UTC (rev 55794) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-18 08:40:42 UTC (rev 55795) @@ -218,7 +218,7 @@ double tolerance) { if (tolerance <= 0.0) - tolerance = 0.01; + tolerance = ON_ZERO_TOLERANCE; if (pointA.DistanceTo(pointB) <= tolerance) { ON_PX_EVENT Event; @@ -234,10 +234,18 @@ } -#define MAX_PCI_DEPTH 8 /** * Point-curve intersections (PCI) */ + +// The maximal depth for subdivision - trade-off between accurancy and +// performance +#define MAX_PCI_DEPTH 8 + +// We make the default tolerance for PSI the same as that of curve-curve +// intersections defined by openNURBS (see other/openNURBS/opennurbs_curve.h) +#define PCI_DEFAULT_TOLERANCE 0.001 + bool ON_Intersect(const ON_3dPoint& pointA, const ON_Curve& curveB, @@ -246,7 +254,7 @@ const ON_Interval* curveB_domain) { if (tolerance <= 0.0) - tolerance = 0.01; + tolerance = PCI_DEFAULT_TOLERANCE; Subcurve root; if (curveB_domain == NULL) { @@ -332,6 +340,11 @@ /** * Point-surface intersections (PSI) */ + +// We make the default tolerance for PSI the same as that of curve-surface +// intersections defined by openNURBS (see other/openNURBS/opennurbs_curve.h) +#define PSI_DEFAULT_TOLERANCE 0.001 + bool ON_Intersect(const ON_3dPoint& pointA, const ON_Surface& surfaceB, @@ -341,7 +354,7 @@ const ON_Interval* surfaceB_vdomain) { if (tolerance <= 0.0) - tolerance = 0.01; + tolerance = PSI_DEFAULT_TOLERANCE; ON_Interval u_domain, v_domain; if (surfaceB_udomain == 0) @@ -542,7 +555,10 @@ }; +// The maximal depth for subdivision - trade-off between accurancy and +// performance #define MAX_SSI_DEPTH 8 + int ON_Intersect(const ON_Surface* surfA, const ON_Surface* surfB, This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-18 08:49:00
|
Revision: 55796 http://sourceforge.net/p/brlcad/code/55796 Author: phoenixyjll Date: 2013-06-18 08:48:54 +0000 (Tue, 18 Jun 2013) Log Message: ----------- The input u_domain and v_domain should be considered. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-18 08:40:42 UTC (rev 55795) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-18 08:48:54 UTC (rev 55796) @@ -257,7 +257,7 @@ tolerance = PCI_DEFAULT_TOLERANCE; Subcurve root; - if (curveB_domain == NULL) { + if (curveB_domain == NULL || *curveB_domain == curveB.Domain()) { root.m_curve = curveB.Duplicate(); root.m_t = curveB.Domain(); } @@ -378,7 +378,9 @@ ON_3dPoint closest_point_3d = surfaceB.PointAt(closest_point_uv.x, closest_point_uv.y); bu_log("%lf %lf %lf\n", pointA.x, pointA.y, pointA.z); bu_log("%lf %lf %lf\n", closest_point_3d.x, closest_point_3d.y, closest_point_3d.z); - if (pointA.DistanceTo(closest_point_3d) <= tolerance) { + if (pointA.DistanceTo(closest_point_3d) <= tolerance + && u_domain.Includes(closest_point_uv.x) + && v_domain.Includes(closest_point_uv.y)) { ON_PX_EVENT Event; Event.m_type = ON_PX_EVENT::psx_point; Event.m_A = pointA; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-20 12:57:59
|
Revision: 55806 http://sourceforge.net/p/brlcad/code/55806 Author: phoenixyjll Date: 2013-06-20 12:57:54 +0000 (Thu, 20 Jun 2013) Log Message: ----------- Use a smaller depth for PSI to improve performance. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-20 12:37:29 UTC (rev 55805) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-20 12:57:54 UTC (rev 55806) @@ -345,6 +345,12 @@ // intersections defined by openNURBS (see other/openNURBS/opennurbs_curve.h) #define PSI_DEFAULT_TOLERANCE 0.001 +// The default maximal depth for creating a surfacee tree (8) is way too +// much - killing performance. Since it's only used for getting an +// estimation, and we use Newton iterations afterwards, it's reasonable +// to use a smaller depth in this step. +#define MAX_PSI_DEPTH 4 + bool ON_Intersect(const ON_3dPoint& pointA, const ON_Surface& surfaceB, @@ -372,16 +378,17 @@ brep->AddSurface(surfaceB.Duplicate()); brep->NewFace(0); ON_2dPoint closest_point_uv; - if (brlcad::get_closest_point(closest_point_uv, brep->Face(0), pointA) == false) { + brlcad::SurfaceTree *tree = new brlcad::SurfaceTree(brep->Face(0), true, MAX_PSI_DEPTH); + if (brlcad::get_closest_point(closest_point_uv, brep->Face(0), pointA, tree) == false) { delete brep; + delete tree; return false; } delete brep; + delete tree; ON_3dPoint closest_point_3d = surfaceB.PointAt(closest_point_uv.x, closest_point_uv.y); - bu_log("%lf %lf %lf\n", pointA.x, pointA.y, pointA.z); - bu_log("%lf %lf %lf\n", closest_point_3d.x, closest_point_3d.y, closest_point_3d.z); if (pointA.DistanceTo(closest_point_3d) <= tolerance && u_domain.Includes(closest_point_uv.x) && v_domain.Includes(closest_point_uv.y)) { This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-21 10:37:29
|
Revision: 55810 http://sourceforge.net/p/brlcad/code/55810 Author: phoenixyjll Date: 2013-06-21 10:37:25 +0000 (Fri, 21 Jun 2013) Log Message: ----------- Add Newton-Raphson iteration in PCI to improve accuracy (use the one from subdivision and linear approximation as a starting point) Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-21 00:09:55 UTC (rev 55809) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-21 10:37:25 UTC (rev 55810) @@ -246,6 +246,11 @@ // intersections defined by openNURBS (see other/openNURBS/opennurbs_curve.h) #define PCI_DEFAULT_TOLERANCE 0.001 +// Newton-Raphson iteration needs a good start. Although we can almost get +// a good starting point every time, but in case of some unfortunate cases, +// we define a maximum iteration, preventing from infinite loop. +#define PCI_MAX_ITERATIONS 100 + bool ON_Intersect(const ON_3dPoint& pointA, const ON_Curve& curveB, @@ -305,29 +310,60 @@ } for (unsigned int i = 0; i < candidates.size(); i++) { - // Use linear approximation to get the intersection point + // Use linear approximation to get an estimated intersection point ON_Line line(candidates[i]->m_curve->PointAtStart(), candidates[i]->m_curve->PointAtEnd()); double t, dis; line.ClosestPointTo(pointA, &t); - // Get the closest point on the line *segment* to the point - double closest_point_t; + + // Make sure line_t belongs to [0, 1] + double line_t; if (t < 0) - closest_point_t = 0; + line_t = 0; else if (t > 1) - closest_point_t = 1; + line_t = 1; else - closest_point_t = t; - ON_3dPoint closest_point = line.PointAt(closest_point_t); - dis = pointA.DistanceTo(closest_point); + line_t = t; + double closest_point_t = candidates[i]->m_t.Min() + candidates[i]->m_t.Length()*line_t; + ON_3dPoint closest_point = curveB.PointAt(closest_point_t); + dis = closest_point.DistanceTo(pointA); + + // Use Newton-Raphson method with the starting point from linear + // approximation. + // It's similar to point projection (or point inversion, or get + // closest point on curve). We calculate the root of: + // (curve(t) - pointA) \dot tangent(t) = 0 (*) + // Then curve(t) may be the closest point on the curve to pointA. + // + // In some cases, the closest point (intersection point) is the + // endpoints of the curve, it does not satisfy (*), but that's + // OK because it already comes out with the linear approximation + // above. + + double last_t = DBL_MAX; + int iterations = 0; + while (!ON_NearZero(last_t-closest_point_t) + && dis > tolerance + && iterations++ < PCI_MAX_ITERATIONS) { + ON_3dVector tangent, s_deriv; + last_t = closest_point_t; + curveB.Ev2Der(closest_point_t, closest_point, tangent, s_deriv); + double F = ON_DotProduct(closest_point-pointA, tangent); + double derivF = ON_DotProduct(closest_point-pointA, s_deriv) + ON_DotProduct(tangent, tangent); + if (!ON_NearZero(derivF)) + closest_point_t -= F/derivF; + dis = closest_point.DistanceTo(pointA); + } + closest_point = curveB.PointAt(closest_point_t); + if (dis <= tolerance) { ON_PX_EVENT *Event = new ON_PX_EVENT; Event->m_type = ON_PX_EVENT::pcx_point; Event->m_A = pointA; Event->m_B = closest_point; - Event->m_b[0] = candidates[i]->m_t.Min()+candidates[i]->m_t.Length()*closest_point_t; + Event->m_b[0] = closest_point_t; Event->m_Mid = (pointA + Event->m_B) * 0.5; - Event->m_radius = dis * 0.5; + Event->m_radius = closest_point.DistanceTo(pointA) * 0.5; x.Append(*Event); return true; } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-24 15:38:42
|
Revision: 55822 http://sourceforge.net/p/brlcad/code/55822 Author: phoenixyjll Date: 2013-06-24 15:38:39 +0000 (Mon, 24 Jun 2013) Log Message: ----------- We have to set removeTrimmed to false, otherwise the surface tree cannot be built correctly. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-24 13:13:46 UTC (rev 55821) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-24 15:38:39 UTC (rev 55822) @@ -414,7 +414,7 @@ brep->AddSurface(surfaceB.Duplicate()); brep->NewFace(0); ON_2dPoint closest_point_uv; - brlcad::SurfaceTree *tree = new brlcad::SurfaceTree(brep->Face(0), true, MAX_PSI_DEPTH); + brlcad::SurfaceTree *tree = new brlcad::SurfaceTree(brep->Face(0), false, MAX_PSI_DEPTH); if (brlcad::get_closest_point(closest_point_uv, brep->Face(0), pointA, tree) == false) { delete brep; delete tree; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-25 07:32:31
|
Revision: 55834 http://sourceforge.net/p/brlcad/code/55834 Author: phoenixyjll Date: 2013-06-25 07:32:26 +0000 (Tue, 25 Jun 2013) Log Message: ----------- overlap_tolerance is used now. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 06:48:22 UTC (rev 55833) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 07:32:26 UTC (rev 55834) @@ -542,7 +542,7 @@ const ON_Curve* curveB, ON_SimpleArray<ON_X_EVENT>& x, double intersection_tolerance, - double UNUSED(overlap_tolerance), + double overlap_tolerance, const ON_Interval* curveA_domain, const ON_Interval* curveB_domain) { This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-25 08:54:28
|
Revision: 55836 http://sourceforge.net/p/brlcad/code/55836 Author: phoenixyjll Date: 2013-06-25 08:54:23 +0000 (Tue, 25 Jun 2013) Log Message: ----------- Fix a fatal typo. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 07:43:03 UTC (rev 55835) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 08:54:23 UTC (rev 55836) @@ -576,10 +576,10 @@ splittedA.push_back((*i).first->m_children[1]); } if ((*i).second->m_islinear || (*i).second->Split() != 0) { - splittedA.push_back((*i).second); + splittedB.push_back((*i).second); } else { - splittedA.push_back((*i).second->m_children[0]); - splittedA.push_back((*i).second->m_children[1]); + splittedB.push_back((*i).second->m_children[0]); + splittedB.push_back((*i).second->m_children[1]); } for (unsigned int j = 0; j < splittedA.size(); j++) for (unsigned int k = 0; k < splittedB.size(); k++) @@ -643,6 +643,7 @@ ON_3dPoint pointB2 = curveB->PointAt(t_b2); double distance1 = pointA1.DistanceTo(pointB1); double distance2 = pointA2.DistanceTo(pointB2); + // Check the validity of the solution if (distance1 < intersection_tolerance && distance2 < intersection_tolerance) { ON_X_EVENT Event; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-25 09:00:30
|
Revision: 55837 http://sourceforge.net/p/brlcad/code/55837 Author: phoenixyjll Date: 2013-06-25 09:00:25 +0000 (Tue, 25 Jun 2013) Log Message: ----------- Add logic in case that one side may fail. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 08:54:23 UTC (rev 55836) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 09:00:25 UTC (rev 55837) @@ -669,6 +669,24 @@ else Event.m_type = ON_X_EVENT::ccx_overlap; x.Append(Event); + } else if (distance1 < intersection_tolerance) { + // in case that the second one was not correct + ON_X_EVENT Event; + Event.m_A[0] = pointA1; + Event.m_B[0] = pointB1; + Event.m_a[0] = t_a1; + Event.m_b[0] = t_b1; + Event.m_type = ON_X_EVENT::ccx_point; + x.Append(Event); + } else if (distance2 < intersection_tolerance) { + // in case that the first one was not correct + ON_X_EVENT Event; + Event.m_A[0] = pointA2; + Event.m_B[0] = pointB2; + Event.m_a[0] = t_a2; + Event.m_b[0] = t_b2; + Event.m_type = ON_X_EVENT::ccx_point; + x.Append(Event); } } } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-25 09:40:18
|
Revision: 55838 http://sourceforge.net/p/brlcad/code/55838 Author: phoenixyjll Date: 2013-06-25 09:40:14 +0000 (Tue, 25 Jun 2013) Log Message: ----------- Fix an opposite logic, and use pointers for ON_X_EVENTs. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 09:00:25 UTC (rev 55837) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 09:40:14 UTC (rev 55838) @@ -524,7 +524,7 @@ J[1][1] = -derivB.y; F[0][0] = pointA.x - pointB.x; F[1][0] = pointA.y - pointB.y; - if (J.Invert(0.0)) { + if (!J.Invert(0.0)) { // FIXME: More elegant error handling. bu_log("Inverse failed.\n"); continue; @@ -626,13 +626,13 @@ double distance = pointA.DistanceTo(pointB); // Check the validity of the solution if (distance < intersection_tolerance) { - ON_X_EVENT Event; - Event.m_A[0] = pointA; - Event.m_B[0] = pointB; - Event.m_a[0] = t_a1; - Event.m_b[0] = t_b1; - Event.m_type = ON_X_EVENT::ccx_point; - x.Append(Event); + ON_X_EVENT *Event = new ON_X_EVENT; + Event->m_A[0] = pointA; + Event->m_B[0] = pointB; + Event->m_a[0] = t_a1; + Event->m_b[0] = t_b1; + Event->m_type = ON_X_EVENT::ccx_point; + x.Append(*Event); } } else { // Check overlap @@ -646,15 +646,15 @@ // Check the validity of the solution if (distance1 < intersection_tolerance && distance2 < intersection_tolerance) { - ON_X_EVENT Event; - Event.m_A[0] = pointA1; - Event.m_A[1] = pointA2; - Event.m_B[0] = pointB1; - Event.m_B[1] = pointB2; - Event.m_a[0] = t_a1; - Event.m_a[1] = t_a2; - Event.m_b[0] = t_b1; - Event.m_b[1] = t_b2; + ON_X_EVENT *Event = new ON_X_EVENT; + Event->m_A[0] = pointA1; + Event->m_A[1] = pointA2; + Event->m_B[0] = pointB1; + Event->m_B[1] = pointB2; + Event->m_a[0] = t_a1; + Event->m_a[1] = t_a2; + Event->m_b[0] = t_b1; + Event->m_b[1] = t_b2; int j; for (j = 1; j < CCI_OVERLAP_TEST_POINTS; j++) { double strike = 1.0/CCI_OVERLAP_TEST_POINTS; @@ -665,28 +665,28 @@ break; } if (j != CCI_OVERLAP_TEST_POINTS) - Event.m_type = ON_X_EVENT::ccx_point; + Event->m_type = ON_X_EVENT::ccx_point; else - Event.m_type = ON_X_EVENT::ccx_overlap; - x.Append(Event); + Event->m_type = ON_X_EVENT::ccx_overlap; + x.Append(*Event); } else if (distance1 < intersection_tolerance) { // in case that the second one was not correct - ON_X_EVENT Event; - Event.m_A[0] = pointA1; - Event.m_B[0] = pointB1; - Event.m_a[0] = t_a1; - Event.m_b[0] = t_b1; - Event.m_type = ON_X_EVENT::ccx_point; - x.Append(Event); + ON_X_EVENT *Event = new ON_X_EVENT; + Event->m_A[0] = pointA1; + Event->m_B[0] = pointB1; + Event->m_a[0] = t_a1; + Event->m_b[0] = t_b1; + Event->m_type = ON_X_EVENT::ccx_point; + x.Append(*Event); } else if (distance2 < intersection_tolerance) { // in case that the first one was not correct - ON_X_EVENT Event; - Event.m_A[0] = pointA2; - Event.m_B[0] = pointB2; - Event.m_a[0] = t_a2; - Event.m_b[0] = t_b2; - Event.m_type = ON_X_EVENT::ccx_point; - x.Append(Event); + ON_X_EVENT *Event = new ON_X_EVENT; + Event->m_A[0] = pointA2; + Event->m_B[0] = pointB2; + Event->m_a[0] = t_a2; + Event->m_b[0] = t_b2; + Event->m_type = ON_X_EVENT::ccx_point; + x.Append(*Event); } } } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-26 02:59:13
|
Revision: 55847 http://sourceforge.net/p/brlcad/code/55847 Author: phoenixyjll Date: 2013-06-26 02:59:09 +0000 (Wed, 26 Jun 2013) Log Message: ----------- Add tolerance in the bounding box intersections. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-25 21:52:09 UTC (rev 55846) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 02:59:09 UTC (rev 55847) @@ -102,10 +102,12 @@ ON_BoundingBox new_bbox(m_node.m_min-vtol, m_node.m_max+vtol); return new_bbox.IsPointIn(pt); } - bool Intersect(const Subcurve& other) const + bool Intersect(const Subcurve& other, double tolerance = 0.0) const { + ON_3dVector vtol(tolerance,tolerance,tolerance); + ON_BoundingBox new_bbox(m_node.m_min-vtol, m_node.m_max+vtol); ON_BoundingBox intersection; - return intersection.Intersection(m_node, other.m_node); + return intersection.Intersection(new_bbox, other.m_node); } }; @@ -559,7 +561,7 @@ typedef std::vector<std::pair<Subcurve*, Subcurve*> > NodePairs; NodePairs candidates, next_candidates; - if (rootA.Intersect(rootB)) + if (rootA.Intersect(rootB, intersection_tolerance)) candidates.push_back(std::make_pair(&rootA, &rootB)); // Use sub-division and bounding box intersections first. @@ -583,7 +585,7 @@ } for (unsigned int j = 0; j < splittedA.size(); j++) for (unsigned int k = 0; k < splittedB.size(); k++) - if (splittedA[j]->Intersect(*splittedB[k])) + if (splittedA[j]->Intersect(*splittedB[k], intersection_tolerance)) next_candidates.push_back(std::make_pair(splittedA[j], splittedB[k])); } candidates = next_candidates; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-26 03:10:26
|
Revision: 55848 http://sourceforge.net/p/brlcad/code/55848 Author: phoenixyjll Date: 2013-06-26 03:10:22 +0000 (Wed, 26 Jun 2013) Log Message: ----------- Check duplication before appending to the array x. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 02:59:09 UTC (rev 55847) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 03:10:22 UTC (rev 55848) @@ -591,6 +591,7 @@ candidates = next_candidates; } + ON_SimpleArray<ON_X_EVENT> tmp_x; // For intersected bounding boxes, we calculate an accurate intersection // point. for (NodePairs::iterator i = candidates.begin(); i != candidates.end(); i++) { @@ -634,7 +635,7 @@ Event->m_a[0] = t_a1; Event->m_b[0] = t_b1; Event->m_type = ON_X_EVENT::ccx_point; - x.Append(*Event); + tmp_x.Append(*Event); } } else { // Check overlap @@ -670,7 +671,7 @@ Event->m_type = ON_X_EVENT::ccx_point; else Event->m_type = ON_X_EVENT::ccx_overlap; - x.Append(*Event); + tmp_x.Append(*Event); } else if (distance1 < intersection_tolerance) { // in case that the second one was not correct ON_X_EVENT *Event = new ON_X_EVENT; @@ -679,7 +680,7 @@ Event->m_a[0] = t_a1; Event->m_b[0] = t_b1; Event->m_type = ON_X_EVENT::ccx_point; - x.Append(*Event); + tmp_x.Append(*Event); } else if (distance2 < intersection_tolerance) { // in case that the first one was not correct ON_X_EVENT *Event = new ON_X_EVENT; @@ -688,11 +689,21 @@ Event->m_a[0] = t_a2; Event->m_b[0] = t_b2; Event->m_type = ON_X_EVENT::ccx_point; - x.Append(*Event); + tmp_x.Append(*Event); } } } + for (int i = 0; i < tmp_x.Count(); i++) { + int j; + for (j = 0; j < x.Count(); j++) { + if (ON_NearZero(tmp_x[i].m_a[0] - x[j].m_a[0]) && ON_NearZero(tmp_x[i].m_a[1] - x[j].m_a[1]) + && ON_NearZero(tmp_x[i].m_b[0] - x[j].m_b[0]) && ON_NearZero(tmp_x[i].m_b[1] - x[j].m_b[1])) + break; + } + if (j == x.Count()) + x.Append(tmp_x[i]); + } return x.Count(); } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-26 04:22:09
|
Revision: 55849 http://sourceforge.net/p/brlcad/code/55849 Author: phoenixyjll Date: 2013-06-26 04:22:06 +0000 (Wed, 26 Jun 2013) Log Message: ----------- If the inverse fails, we try another two directions. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 03:10:22 UTC (rev 55848) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 04:22:06 UTC (rev 55849) @@ -505,7 +505,7 @@ // y_a(t_a) - y_b(t_b) = 0 // z_a(t_a) - z_b(t_b) = 0 // It's an over-determined system. - // We use Newton-Raphson iterations to solve the first two equations, + // We use Newton-Raphson iterations to solve two equations of the three, // and use the third for checking. double last_t_a = DBL_MAX*.5, last_t_b = DBL_MAX*.5; ON_3dPoint pointA = curveA->PointAt(t_a); @@ -527,9 +527,26 @@ F[0][0] = pointA.x - pointB.x; F[1][0] = pointA.y - pointB.y; if (!J.Invert(0.0)) { - // FIXME: More elegant error handling. - bu_log("Inverse failed.\n"); - continue; + // bu_log("Inverse failed.\n"); + J[0][0] = derivA.x; + J[0][1] = -derivB.x; + J[1][0] = derivA.z; + J[1][1] = -derivB.z; + F[0][0] = pointA.x - pointB.x; + F[1][0] = pointA.z - pointB.z; + if (!J.Invert(0.0)) { + // bu_log("Inverse failed again.\n"); + J[0][0] = derivA.y; + J[0][1] = -derivB.y; + J[1][0] = derivA.z; + J[1][1] = -derivB.z; + F[0][0] = pointA.y - pointB.y; + F[1][0] = pointA.z - pointB.z; + if (!J.Invert(0.0)) { + // bu_log("Inverse failed and never try again.\n"); + continue; + } + } } ON_Matrix Delta; Delta.Multiply(J, F); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-26 09:17:31
|
Revision: 55851 http://sourceforge.net/p/brlcad/code/55851 Author: phoenixyjll Date: 2013-06-26 09:17:26 +0000 (Wed, 26 Jun 2013) Log Message: ----------- More work on the tolerance value to get a more accurate and correct result. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 09:09:46 UTC (rev 55850) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 09:17:26 UTC (rev 55851) @@ -498,7 +498,7 @@ } void -newton_cci(double& t_a, double& t_b, const ON_Curve* curveA, const ON_Curve* curveB, double intersection_tolerance) +newton_cci(double& t_a, double& t_b, const ON_Curve* curveA, const ON_Curve* curveB) { // Equations: // x_a(t_a) - x_b(t_b) = 0 @@ -513,7 +513,6 @@ double distance = pointA.DistanceTo(pointB); int iteration = 0; while (fabs(last_t_a - t_a) + fabs(last_t_b - t_b) > ON_ZERO_TOLERANCE - && distance > intersection_tolerance && iteration++ < CCI_MAX_ITERATIONS) { last_t_a = t_a, last_t_b = t_b; ON_3dVector derivA, derivB; @@ -634,12 +633,16 @@ // If they converge to one point, it's considered an intersection // point, otherwise it's considered an overlap event. double t_a1 = i->first->m_t.Min(), t_b1 = i->second->m_t.Min(); - newton_cci(t_a1, t_b1, curveA, curveB, intersection_tolerance); + newton_cci(t_a1, t_b1, curveA, curveB); double t_a2 = i->first->m_t.Max(), t_b2 = i->second->m_t.Max(); - newton_cci(t_a2, t_b2, curveA, curveB, intersection_tolerance); + newton_cci(t_a2, t_b2, curveA, curveB); - if (fabs(t_a1 - t_a2) + fabs(t_b1 - t_b2) < ON_ZERO_TOLERANCE*2) { - // the errors may double, so we use ON_ZERO_TOLERANCE*2 + ON_3dPoint pointA1 = curveA->PointAt(t_a1); + ON_3dPoint pointB1 = curveB->PointAt(t_b1); + ON_3dPoint pointA2 = curveA->PointAt(t_a2); + ON_3dPoint pointB2 = curveB->PointAt(t_b2); + if (pointA1.DistanceTo(pointA2) < intersection_tolerance + && pointB1.DistanceTo(pointB2) < intersection_tolerance) { // it's considered the same point ON_3dPoint pointA = curveA->PointAt(t_a1); ON_3dPoint pointB = curveB->PointAt(t_b1); @@ -656,11 +659,7 @@ } } else { // Check overlap - bu_log("Maybe overlap.\n"); - ON_3dPoint pointA1 = curveA->PointAt(t_a1); - ON_3dPoint pointB1 = curveB->PointAt(t_b1); - ON_3dPoint pointA2 = curveA->PointAt(t_a2); - ON_3dPoint pointB2 = curveB->PointAt(t_b2); + // bu_log("Maybe overlap.\n"); double distance1 = pointA1.DistanceTo(pointB1); double distance2 = pointA2.DistanceTo(pointB2); @@ -714,8 +713,10 @@ for (int i = 0; i < tmp_x.Count(); i++) { int j; for (j = 0; j < x.Count(); j++) { - if (ON_NearZero(tmp_x[i].m_a[0] - x[j].m_a[0]) && ON_NearZero(tmp_x[i].m_a[1] - x[j].m_a[1]) - && ON_NearZero(tmp_x[i].m_b[0] - x[j].m_b[0]) && ON_NearZero(tmp_x[i].m_b[1] - x[j].m_b[1])) + if (tmp_x[i].m_A[0].DistanceTo(x[j].m_A[0]) < intersection_tolerance + && tmp_x[i].m_A[1].DistanceTo(x[j].m_A[1]) < intersection_tolerance + && tmp_x[i].m_B[0].DistanceTo(x[j].m_B[0]) < intersection_tolerance + && tmp_x[i].m_B[1].DistanceTo(x[j].m_B[1]) < intersection_tolerance) break; } if (j == x.Count()) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-26 09:21:05
|
Revision: 55852 http://sourceforge.net/p/brlcad/code/55852 Author: phoenixyjll Date: 2013-06-26 09:21:02 +0000 (Wed, 26 Jun 2013) Log Message: ----------- Remove the set-but-not-used variable distance. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 09:17:26 UTC (rev 55851) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 09:21:02 UTC (rev 55852) @@ -510,7 +510,6 @@ double last_t_a = DBL_MAX*.5, last_t_b = DBL_MAX*.5; ON_3dPoint pointA = curveA->PointAt(t_a); ON_3dPoint pointB = curveB->PointAt(t_b); - double distance = pointA.DistanceTo(pointB); int iteration = 0; while (fabs(last_t_a - t_a) + fabs(last_t_b - t_b) > ON_ZERO_TOLERANCE && iteration++ < CCI_MAX_ITERATIONS) { @@ -551,7 +550,6 @@ Delta.Multiply(J, F); t_a -= Delta[0][0]; t_b -= Delta[1][0]; - distance = pointA.DistanceTo(pointB); } } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <br...@us...> - 2013-06-26 16:01:40
|
Revision: 55857 http://sourceforge.net/p/brlcad/code/55857 Author: brlcad Date: 2013-06-26 16:01:37 +0000 (Wed, 26 Jun 2013) Log Message: ----------- style conformance Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 15:58:30 UTC (rev 55856) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 16:01:37 UTC (rev 55857) @@ -98,13 +98,13 @@ } bool IsPointIn(const ON_3dPoint &pt, double tolerance = 0.0) { - ON_3dVector vtol(tolerance,tolerance,tolerance); + ON_3dVector vtol(tolerance, tolerance, tolerance); ON_BoundingBox new_bbox(m_node.m_min-vtol, m_node.m_max+vtol); return new_bbox.IsPointIn(pt); } bool Intersect(const Subcurve& other, double tolerance = 0.0) const { - ON_3dVector vtol(tolerance,tolerance,tolerance); + ON_3dVector vtol(tolerance, tolerance, tolerance); ON_BoundingBox new_bbox(m_node.m_min-vtol, m_node.m_max+vtol); ON_BoundingBox intersection; return intersection.Intersection(new_bbox, other.m_node); @@ -497,6 +497,7 @@ return true; } + void newton_cci(double& t_a, double& t_b, const ON_Curve* curveA, const ON_Curve* curveB) { @@ -553,6 +554,7 @@ } } + int ON_Intersect(const ON_Curve* curveA, const ON_Curve* curveB, @@ -963,18 +965,22 @@ if ((*i).second->m_children[0] == NULL) ret2 = (*i).second->Split(); if (ret1) { - if (ret2) { /* both splits failed */ + if (ret2) { + /* both splits failed */ tmp_pairs.push_back(*i); h = MAX_SSI_DEPTH; - } else { /* the first failed */ + } else { + /* the first failed */ for (int j = 0; j < 4; j++) tmp_pairs.push_back(std::make_pair((*i).first, (*i).second->m_children[j])); } } else { - if (ret2) { /* the second failed */ + if (ret2) { + /* the second failed */ for (int j = 0; j < 4; j++) tmp_pairs.push_back(std::make_pair((*i).first->m_children[j], (*i).second)); - } else { /* both success */ + } else { + /* both success */ for (int j = 0; j < 4; j++) for (int k = 0; k < 4; k++) tmp_pairs.push_back(std::make_pair((*i).first->m_children[j], (*i).second->m_children[k])); @@ -1294,6 +1300,7 @@ return x.Count(); } + // Local Variables: // tab-width: 8 // mode: C++ This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-27 05:32:50
|
Revision: 55866 http://sourceforge.net/p/brlcad/code/55866 Author: phoenixyjll Date: 2013-06-27 05:32:46 +0000 (Thu, 27 Jun 2013) Log Message: ----------- Merge the overlap events that are continuous, and eliminate the intersection points that are inside the overlap events. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-26 19:07:29 UTC (rev 55865) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-27 05:32:46 UTC (rev 55866) @@ -552,10 +552,23 @@ t_a -= Delta[0][0]; t_b -= Delta[1][0]; } + + // Make sure the solution is inside the domains + t_a = std::min(t_a, curveA->Domain().Max()); + t_a = std::max(t_a, curveA->Domain().Min()); + t_b = std::min(t_b, curveB->Domain().Max()); + t_b = std::max(t_b, curveB->Domain().Min()); } int +compare_by_m_a0(const ON_X_EVENT* a, const ON_X_EVENT* b) +{ + return a->m_a[0] - b->m_a[0]; +} + + +int ON_Intersect(const ON_Curve* curveA, const ON_Curve* curveB, ON_SimpleArray<ON_X_EVENT>& x, @@ -666,14 +679,26 @@ // Check the validity of the solution if (distance1 < intersection_tolerance && distance2 < intersection_tolerance) { ON_X_EVENT *Event = new ON_X_EVENT; - Event->m_A[0] = pointA1; - Event->m_A[1] = pointA2; - Event->m_B[0] = pointB1; - Event->m_B[1] = pointB2; - Event->m_a[0] = t_a1; - Event->m_a[1] = t_a2; - Event->m_b[0] = t_b1; - Event->m_b[1] = t_b2; + // We make sure that m_a[0] <= m_a[1] + if (t_a1 <= t_a2) { + Event->m_A[0] = pointA1; + Event->m_A[1] = pointA2; + Event->m_B[0] = pointB1; + Event->m_B[1] = pointB2; + Event->m_a[0] = t_a1; + Event->m_a[1] = t_a2; + Event->m_b[0] = t_b1; + Event->m_b[1] = t_b2; + } else { + Event->m_A[0] = pointA2; + Event->m_A[1] = pointA1; + Event->m_B[0] = pointB2; + Event->m_B[1] = pointB1; + Event->m_a[0] = t_a2; + Event->m_a[1] = t_a1; + Event->m_b[0] = t_b2; + Event->m_b[1] = t_b1; + } int j; for (j = 1; j < CCI_OVERLAP_TEST_POINTS; j++) { double strike = 1.0/CCI_OVERLAP_TEST_POINTS; @@ -710,18 +735,66 @@ } } + ON_SimpleArray<ON_X_EVENT> points, overlap, pending; + // Use an O(n^2) naive approach to eliminate duplications for (int i = 0; i < tmp_x.Count(); i++) { int j; - for (j = 0; j < x.Count(); j++) { - if (tmp_x[i].m_A[0].DistanceTo(x[j].m_A[0]) < intersection_tolerance - && tmp_x[i].m_A[1].DistanceTo(x[j].m_A[1]) < intersection_tolerance - && tmp_x[i].m_B[0].DistanceTo(x[j].m_B[0]) < intersection_tolerance - && tmp_x[i].m_B[1].DistanceTo(x[j].m_B[1]) < intersection_tolerance) + if (tmp_x[i].m_type == ON_X_EVENT::ccx_overlap) { + overlap.Append(tmp_x[i]); + continue; + } + // ccx_point + for (j = 0; j < points.Count(); j++) { + if (tmp_x[i].m_A[0].DistanceTo(points[j].m_A[0]) < intersection_tolerance + && tmp_x[i].m_A[1].DistanceTo(points[j].m_A[1]) < intersection_tolerance + && tmp_x[i].m_B[0].DistanceTo(points[j].m_B[0]) < intersection_tolerance + && tmp_x[i].m_B[1].DistanceTo(points[j].m_B[1]) < intersection_tolerance) break; } - if (j == x.Count()) - x.Append(tmp_x[i]); + if (j == points.Count()) { + points.Append(tmp_x[i]); + } } + + // Merge the overlap events that are continuous + overlap.QuickSort(compare_by_m_a0); + for (int i = 0; i < overlap.Count(); i++) { + bool merged = false; + for (int j = 0; j < pending.Count(); j++) { + if (pending[j].m_a[1] < overlap[i].m_a[0] - intersection_tolerance) { + x.Append(pending[j]); + pending.Remove(j); + j--; + continue; + } + if (overlap[i].m_a[0] < pending[j].m_a[1] + intersection_tolerance + && overlap[i].m_b[0] < pending[j].m_b[1] + intersection_tolerance) { + // Need to merge + merged = true; + pending[j].m_a[1] = std::max(overlap[i].m_a[1], pending[j].m_a[1]); + pending[j].m_b[1] = std::max(overlap[i].m_b[1], pending[j].m_b[1]); + break; + } + } + if (merged == false) + pending.Append(overlap[i]); + } + for (int i = 0; i < pending.Count(); i++) + x.Append(pending[i]); + + // The intersection points shouldn't be inside the overlapped parts. + int overlap_events = x.Count(); + for (int i = 0; i < points.Count(); i++) { + int j; + for (j = 0; j < overlap_events; j++) { + if (points[i].m_a[0] > x[j].m_a[0] - intersection_tolerance + && points[i].m_a[0] < x[j].m_a[1] + intersection_tolerance) + break; + } + if (j == overlap_events) + x.Append(points[i]); + } + return x.Count(); } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-27 09:30:58
|
Revision: 55869 http://sourceforge.net/p/brlcad/code/55869 Author: phoenixyjll Date: 2013-06-27 09:30:55 +0000 (Thu, 27 Jun 2013) Log Message: ----------- Some special handling for linear curves. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-27 08:35:08 UTC (rev 55868) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-27 09:30:55 UTC (rev 55869) @@ -624,113 +624,160 @@ // For intersected bounding boxes, we calculate an accurate intersection // point. for (NodePairs::iterator i = candidates.begin(); i != candidates.end(); i++) { - /* - // First, use linear approximation to get a starting point - ON_Line lineA(curveA->PointAt(i->first->m_t.Min()), curveA->PointAt(i->first->m_t.Max())); - ON_Line lineB(curveB->PointAt(i->second->m_t.Min()), curveB->PointAt(i->second->m_t.Max())); - double t_lineA, t_lineB; - double t_a, t_b; - if (ON_IntersectLineLine(lineA, lineB, &t_lineA, &t_lineB, ON_ZERO_TOLERANCE, true)) { - // The line segments intersect - t_a = i->first->m_t.ParameterAt(t_lineA); - t_b = i->second->m_t.ParameterAt(t_lineB); - } else { - // Sometimes the approximated line segments do not intersect, - // but the curves DO intersect. We use the mid-points of the - // sub-curves as the starting point. - t_a = i->first->m_t.Mid(); - t_b = i->second->m_t.Mid(); - }*/ + if (i->first->m_islinear && i->second->m_islinear) { + // Both of them are linear, we use ON_IntersectLineLine + ON_Line lineA(curveA->PointAt(i->first->m_t.Min()), curveA->PointAt(i->first->m_t.Max())); + ON_Line lineB(curveB->PointAt(i->second->m_t.Min()), curveB->PointAt(i->second->m_t.Max())); + if (lineA.Direction().IsParallelTo(lineB.Direction())) { + // They are parallel + if (lineA.MinimumDistanceTo(lineB) < intersection_tolerance) { + // we report a ccx_overlap event + double startB_on_A, endB_on_A, startA_on_B, endA_on_B; + lineA.ClosestPointTo(lineB.from, &startB_on_A); + lineA.ClosestPointTo(lineB.to, &endB_on_A); + lineB.ClosestPointTo(lineA.from, &startA_on_B); + lineB.ClosestPointTo(lineA.to, &endA_on_B); - // Use two different start points - the two end-points of the interval - // If they converge to one point, it's considered an intersection - // point, otherwise it's considered an overlap event. - double t_a1 = i->first->m_t.Min(), t_b1 = i->second->m_t.Min(); - newton_cci(t_a1, t_b1, curveA, curveB); - double t_a2 = i->first->m_t.Max(), t_b2 = i->second->m_t.Max(); - newton_cci(t_a2, t_b2, curveA, curveB); + if (startB_on_A > 1+intersection_tolerance || endB_on_A < -intersection_tolerance + || startA_on_B > 1+intersection_tolerance || endA_on_B < -intersection_tolerance) + continue; - ON_3dPoint pointA1 = curveA->PointAt(t_a1); - ON_3dPoint pointB1 = curveB->PointAt(t_b1); - ON_3dPoint pointA2 = curveA->PointAt(t_a2); - ON_3dPoint pointB2 = curveB->PointAt(t_b2); - if (pointA1.DistanceTo(pointA2) < intersection_tolerance - && pointB1.DistanceTo(pointB2) < intersection_tolerance) { - // it's considered the same point - ON_3dPoint pointA = curveA->PointAt(t_a1); - ON_3dPoint pointB = curveB->PointAt(t_b1); - double distance = pointA.DistanceTo(pointB); - // Check the validity of the solution - if (distance < intersection_tolerance) { - ON_X_EVENT *Event = new ON_X_EVENT; - Event->m_A[0] = pointA; - Event->m_B[0] = pointB; - Event->m_a[0] = t_a1; - Event->m_b[0] = t_b1; - Event->m_type = ON_X_EVENT::ccx_point; - tmp_x.Append(*Event); + double t_a1, t_a2, t_b1, t_b2; + if (startB_on_A > 0.0) + t_a1 = startB_on_A, t_b1 = 0.0; + else + t_a1 = 0.0, t_b1 = startA_on_B; + if (endB_on_A < 1.0) + t_a2 = endB_on_A, t_b2 = 1.0; + else + t_a2 = 1.0, t_b2 = endA_on_B; + + ON_X_EVENT* Event = new ON_X_EVENT; + Event->m_A[0] = lineA.PointAt(t_a1); + Event->m_A[1] = lineA.PointAt(t_a2); + Event->m_B[0] = lineB.PointAt(t_b1); + Event->m_B[1] = lineB.PointAt(t_b2); + Event->m_a[0] = i->first->m_t.ParameterAt(t_a1); + Event->m_a[1] = i->first->m_t.ParameterAt(t_a2); + Event->m_b[0] = i->second->m_t.ParameterAt(t_b1); + Event->m_b[1] = i->second->m_t.ParameterAt(t_b2); + Event->m_type = ON_X_EVENT::TYPE::ccx_overlap; + tmp_x.Append(*Event); + } + } else { + // They are not parallel, check intersection point + double t_lineA, t_lineB; + double t_a, t_b; + if (ON_IntersectLineLine(lineA, lineB, &t_lineA, &t_lineB, ON_ZERO_TOLERANCE, true)) { + // The line segments intersect + t_a = i->first->m_t.ParameterAt(t_lineA); + t_b = i->second->m_t.ParameterAt(t_lineB); + + ON_X_EVENT* Event = new ON_X_EVENT; + Event->m_A[0] = lineA.PointAt(t_lineA); + Event->m_B[0] = lineB.PointAt(t_lineB); + Event->m_a[0] = t_a; + Event->m_b[0] = t_b; + Event->m_type = ON_X_EVENT::TYPE::ccx_point; + tmp_x.Append(*Event); + } } } else { - // Check overlap - // bu_log("Maybe overlap.\n"); - double distance1 = pointA1.DistanceTo(pointB1); - double distance2 = pointA2.DistanceTo(pointB2); + // They are not both linear. - // Check the validity of the solution - if (distance1 < intersection_tolerance && distance2 < intersection_tolerance) { - ON_X_EVENT *Event = new ON_X_EVENT; - // We make sure that m_a[0] <= m_a[1] - if (t_a1 <= t_a2) { + // Use two different start points - the two end-points of the interval + // If they converge to one point, it's considered an intersection + // point, otherwise it's considered an overlap event. + // FIXME: Find a better machanism to check overlapping, because this method + // may miss some overlap cases. (Overlap events can also converge to one + // point) + double t_a1 = i->first->m_t.Min(), t_b1 = i->second->m_t.Min(); + newton_cci(t_a1, t_b1, curveA, curveB); + double t_a2 = i->first->m_t.Max(), t_b2 = i->second->m_t.Max(); + newton_cci(t_a2, t_b2, curveA, curveB); + + ON_3dPoint pointA1 = curveA->PointAt(t_a1); + ON_3dPoint pointB1 = curveB->PointAt(t_b1); + ON_3dPoint pointA2 = curveA->PointAt(t_a2); + ON_3dPoint pointB2 = curveB->PointAt(t_b2); + if (pointA1.DistanceTo(pointA2) < intersection_tolerance + && pointB1.DistanceTo(pointB2) < intersection_tolerance) { + // it's considered the same point + ON_3dPoint pointA = curveA->PointAt(t_a1); + ON_3dPoint pointB = curveB->PointAt(t_b1); + double distance = pointA.DistanceTo(pointB); + // Check the validity of the solution + if (distance < intersection_tolerance) { + ON_X_EVENT *Event = new ON_X_EVENT; + Event->m_A[0] = pointA; + Event->m_B[0] = pointB; + Event->m_a[0] = t_a1; + Event->m_b[0] = t_b1; + Event->m_type = ON_X_EVENT::ccx_point; + tmp_x.Append(*Event); + } + } else { + // Check overlap + // bu_log("Maybe overlap.\n"); + double distance1 = pointA1.DistanceTo(pointB1); + double distance2 = pointA2.DistanceTo(pointB2); + + // Check the validity of the solution + if (distance1 < intersection_tolerance && distance2 < intersection_tolerance) { + ON_X_EVENT *Event = new ON_X_EVENT; + // We make sure that m_a[0] <= m_a[1] + if (t_a1 <= t_a2) { + Event->m_A[0] = pointA1; + Event->m_A[1] = pointA2; + Event->m_B[0] = pointB1; + Event->m_B[1] = pointB2; + Event->m_a[0] = t_a1; + Event->m_a[1] = t_a2; + Event->m_b[0] = t_b1; + Event->m_b[1] = t_b2; + } else { + Event->m_A[0] = pointA2; + Event->m_A[1] = pointA1; + Event->m_B[0] = pointB2; + Event->m_B[1] = pointB1; + Event->m_a[0] = t_a2; + Event->m_a[1] = t_a1; + Event->m_b[0] = t_b2; + Event->m_b[1] = t_b1; + } + int j; + for (j = 1; j < CCI_OVERLAP_TEST_POINTS; j++) { + double strike = 1.0/CCI_OVERLAP_TEST_POINTS; + ON_3dPoint test_point = curveA->PointAt(t_a1 + (t_a2-t_a1)*j*strike); + ON_ClassArray<ON_PX_EVENT> pci_x; + // Use point-curve intersection + if (!ON_Intersect(test_point, *curveA, pci_x, overlap_tolerance)) + break; + } + if (j != CCI_OVERLAP_TEST_POINTS) + Event->m_type = ON_X_EVENT::ccx_point; + else + Event->m_type = ON_X_EVENT::ccx_overlap; + tmp_x.Append(*Event); + } else if (distance1 < intersection_tolerance) { + // in case that the second one was not correct + ON_X_EVENT *Event = new ON_X_EVENT; Event->m_A[0] = pointA1; - Event->m_A[1] = pointA2; Event->m_B[0] = pointB1; - Event->m_B[1] = pointB2; Event->m_a[0] = t_a1; - Event->m_a[1] = t_a2; Event->m_b[0] = t_b1; - Event->m_b[1] = t_b2; - } else { + Event->m_type = ON_X_EVENT::ccx_point; + tmp_x.Append(*Event); + } else if (distance2 < intersection_tolerance) { + // in case that the first one was not correct + ON_X_EVENT *Event = new ON_X_EVENT; Event->m_A[0] = pointA2; - Event->m_A[1] = pointA1; Event->m_B[0] = pointB2; - Event->m_B[1] = pointB1; Event->m_a[0] = t_a2; - Event->m_a[1] = t_a1; Event->m_b[0] = t_b2; - Event->m_b[1] = t_b1; + Event->m_type = ON_X_EVENT::ccx_point; + tmp_x.Append(*Event); } - int j; - for (j = 1; j < CCI_OVERLAP_TEST_POINTS; j++) { - double strike = 1.0/CCI_OVERLAP_TEST_POINTS; - ON_3dPoint test_point = curveA->PointAt(t_a1 + (t_a2-t_a1)*j*strike); - ON_ClassArray<ON_PX_EVENT> pci_x; - // Use point-curve intersection - if (!ON_Intersect(test_point, *curveA, pci_x, overlap_tolerance)) - break; - } - if (j != CCI_OVERLAP_TEST_POINTS) - Event->m_type = ON_X_EVENT::ccx_point; - else - Event->m_type = ON_X_EVENT::ccx_overlap; - tmp_x.Append(*Event); - } else if (distance1 < intersection_tolerance) { - // in case that the second one was not correct - ON_X_EVENT *Event = new ON_X_EVENT; - Event->m_A[0] = pointA1; - Event->m_B[0] = pointB1; - Event->m_a[0] = t_a1; - Event->m_b[0] = t_b1; - Event->m_type = ON_X_EVENT::ccx_point; - tmp_x.Append(*Event); - } else if (distance2 < intersection_tolerance) { - // in case that the first one was not correct - ON_X_EVENT *Event = new ON_X_EVENT; - Event->m_A[0] = pointA2; - Event->m_B[0] = pointB2; - Event->m_a[0] = t_a2; - Event->m_b[0] = t_b2; - Event->m_type = ON_X_EVENT::ccx_point; - tmp_x.Append(*Event); } } } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <pho...@us...> - 2013-06-27 09:39:24
|
Revision: 55870 http://sourceforge.net/p/brlcad/code/55870 Author: phoenixyjll Date: 2013-06-27 09:39:22 +0000 (Thu, 27 Jun 2013) Log Message: ----------- remove trailing ws. Modified Paths: -------------- brlcad/trunk/src/libbrep/intersect.cpp Modified: brlcad/trunk/src/libbrep/intersect.cpp =================================================================== --- brlcad/trunk/src/libbrep/intersect.cpp 2013-06-27 09:30:55 UTC (rev 55869) +++ brlcad/trunk/src/libbrep/intersect.cpp 2013-06-27 09:39:22 UTC (rev 55870) @@ -800,7 +800,7 @@ } if (j == points.Count()) { points.Append(tmp_x[i]); - } + } } // Merge the overlap events that are continuous This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |