Skip to content

Commit

Permalink
fixing warnings and cleanup post-debugging incircle2d() for expansions
Browse files Browse the repository at this point in the history
  • Loading branch information
BrunoLevy committed Dec 8, 2023
1 parent 5c450c7 commit a092052
Show file tree
Hide file tree
Showing 3 changed files with 42 additions and 107 deletions.
78 changes: 38 additions & 40 deletions src/lib/geogram/mesh/mesh_CSG.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -688,8 +688,8 @@ namespace GEO {
z = 1.0 - z;
}
if(center) {
x -= double(nu-1.0)/2.0;
y -= double(nv-1.0)/2.0;
x -= double(nu-1)/2.0;
y -= double(nv-1)/2.0;
}
p[0] = x;
p[1] = y;
Expand Down Expand Up @@ -1435,10 +1435,10 @@ namespace GEO {
C->vertices.point_ptr(v)[0] += (x1 - dx);
C->vertices.point_ptr(v)[1] += (y1 - dy);
}
CSGScope scope;
scope.push_back(result);
scope.push_back(C);
result = difference(scope);
CSGScope scope2;
scope2.push_back(result);
scope2.push_back(C);
result = difference(scope2);
vector<index_t> remove_f(result->facets.nb(),0);
for(index_t f: result->facets) {
for(index_t lv=0; lv<result->facets.nb_vertices(f); ++lv) {
Expand All @@ -1453,48 +1453,46 @@ namespace GEO {
result->facets.compute_borders();
result->vertices.set_dimension(2);
result->update_bbox();
return result;
} else {
// Super unelegant brute force algorithm !!
// But just extracting the silhouette does not work because
// we need a smarter in/out classification algorithm, that
// seemingly cannot be expressed as a boolean expression
// passed to triangulate()
{
CSGScope scope;
for(index_t f: result->facets) {
const double* p1 = result->vertices.point_ptr(
result->facets.vertex(f,0)
);

const double* p2 = result->vertices.point_ptr(
result->facets.vertex(f,1)
);

const double* p3 = result->vertices.point_ptr(
result->facets.vertex(f,2)
);

// I thought that I could have said here: != POSITIVE
// NOTE: different set of isolated vertices each time,
// there is something not normal.
if(PCK::orient_2d(p1,p2,p3) == ZERO) {
continue;
CSGScope scope2;
for(index_t f: result->facets) {
const double* p1 = result->vertices.point_ptr(
result->facets.vertex(f,0)
);

const double* p2 = result->vertices.point_ptr(
result->facets.vertex(f,1)
);

const double* p3 = result->vertices.point_ptr(
result->facets.vertex(f,2)
);

// I thought that I could have said here: != POSITIVE
// NOTE: different set of isolated vertices each time,
// there is something not normal.
if(PCK::orient_2d(p1,p2,p3) == ZERO) {
continue;
}

CSGMesh_var F = new CSGMesh;
F->vertices.set_dimension(2);

F->vertices.create_vertex(p1);
F->vertices.create_vertex(p2);
F->vertices.create_vertex(p3);
F->facets.create_triangle(0,1,2);
F->facets.compute_borders();
F->update_bbox();
scope2.push_back(F);
}

CSGMesh_var F = new CSGMesh;
F->vertices.set_dimension(2);

F->vertices.create_vertex(p1);
F->vertices.create_vertex(p2);
F->vertices.create_vertex(p3);
F->facets.create_triangle(0,1,2);
F->facets.compute_borders();
F->update_bbox();
scope.push_back(F);
}
result = union_instr(scope);
return result;
result = union_instr(scope2);
}
}
return result;
Expand Down
6 changes: 3 additions & 3 deletions src/lib/geogram/mesh/mesh_surface_intersection.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1446,7 +1446,7 @@ namespace GEO {
}
} else {
// Else compute the radial sort geometrically
bool OK = I_.radial_bundles_.radial_sort(bndl, RS);
OK = I_.radial_bundles_.radial_sort(bndl, RS);
// May return !OK (if not using geogram_plus) when it
// cannot sort (example_022.csg and example_024.csg)
if(!OK) {
Expand Down Expand Up @@ -1989,8 +1989,8 @@ namespace GEO {
geo_assert(v1 != index_t(-1));
geo_assert(v2 != index_t(-1));
geo_assert(v3 != index_t(-1));
index_t f = mesh_.facets.create_triangle(v1,v2,v3);
facet_group[f] = current_group;
index_t new_f = mesh_.facets.create_triangle(v1,v2,v3);
facet_group[new_f] = current_group;
}
}
}
Expand Down
65 changes: 1 addition & 64 deletions src/lib/geogram/numerics/exact_geometry.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -328,7 +328,7 @@ namespace GEO {
// worth it.

// Filter
if(false) {
{
interval_nt::Rounding rounding;
interval_nt l3I(l3);
interval_nt L1 = interval_nt(l0) - l3I;
Expand Down Expand Up @@ -372,70 +372,7 @@ namespace GEO {

// Exact
stats.log_exact();

{
rational_nt l0_(
(geo_sqr(p0.x) + geo_sqr(p0.y)).estimate(),
geo_sqr(p0.w).estimate()
);

rational_nt l1_(
(geo_sqr(p1.x) + geo_sqr(p1.y)).estimate(),
geo_sqr(p1.w).estimate()
);

rational_nt l2_(
(geo_sqr(p2.x) + geo_sqr(p2.y)).estimate(),
geo_sqr(p2.w).estimate()
);

rational_nt l3_(
(geo_sqr(p3.x) + geo_sqr(p3.y)).estimate(),
geo_sqr(p3.w).estimate()
);


l0_=rational_nt(l0_.num().estimate()/l0_.denom().estimate(),1.0);
l1_=rational_nt(l1_.num().estimate()/l1_.denom().estimate(),1.0);
l2_=rational_nt(l2_.num().estimate()/l2_.denom().estimate(),1.0);
l3_=rational_nt(l3_.num().estimate()/l3_.denom().estimate(),1.0);

l0_.optimize(); l1_.optimize(); l2_.optimize(); l3_.optimize();


l0_ = rational_nt(l0,1.0);
l1_ = rational_nt(l1,1.0);
l2_ = rational_nt(l2,1.0);
l3_ = rational_nt(l3,1.0);

vec2HE P1 = p0 - p3;
vec2HE P2 = p1 - p3;
vec2HE P3 = p2 - p3;
P1.optimize(); P2.optimize(); P3.optimize();

rational_nt L1 = l0_ - l3_;
rational_nt L2 = l1_ - l3_;
rational_nt L3 = l2_ - l3_;
L1.optimize(); L2.optimize(); L3.optimize();

expansion_nt M1 = det2x2(P2.x, P2.y, P3.x, P3.y);
expansion_nt M2 = det2x2(P1.x, P1.y, P3.x, P3.y);
expansion_nt M3 = det2x2(P1.x, P1.y, P2.x, P2.y);
M1.optimize(); M2.optimize(); M3.optimize();

expansion_nt D =
L1.num() * L2.denom() * L3.denom() * P1.w * M1 -
L2.num() * L1.denom() * L3.denom() * P2.w * M2 +
L3.num() * L1.denom() * L2.denom() * P3.w * M3 ;

result = Sign(
D.sign() *
P1.w.sign() * P2.w.sign() * P3.w.sign() *
L1.denom().sign() * L2.denom().sign() * L3.denom().sign()
);
}

if(false) {
expansion_nt L1(expansion_nt::DIFF, l0, l3);
expansion_nt L2(expansion_nt::DIFF, l1, l3);
expansion_nt L3(expansion_nt::DIFF, l2, l3);
Expand Down

0 comments on commit a092052

Please sign in to comment.