commit 6b6b61c052b077a2c71b70da00ad93f4b4a5357c
parent 6adc14c2906c90094f239c2af67327cd2a2a4b89
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date: Tue, 27 Feb 2024 15:59:29 +0100
BugFix: wrong grouping of connex components
Enclosures are made of CC that are grouped according to topology.
in the first grouping stage, each CC is either grouped to a single other
CC, or this CC is part of enclosure 0 (the external boundary of the
system).
Grouping is made by finding the closest vertex Vc upwards from the
highest vertex Vh of the CC. Vc + the side of the triangle it belongs to
give the CC to group with. If the selected triangle is nearly parallel
to the line of sight, it is not possible to accurately determine which
side is seen, and the grouping can select the wrong CC.
The fix consists in detecting the nearly parallel situations and to add
a postprocess: get all the triangles around Vc and trace a ray from Vh's
barycenter to Vc. The first hit with accurate determination of the hit
side gives the CC.
Diffstat:
1 file changed, 375 insertions(+), 214 deletions(-)
diff --git a/src/senc3d_scene_analyze.c b/src/senc3d_scene_analyze.c
@@ -27,6 +27,7 @@
#include <rsys/mem_allocator.h>
#include <rsys/hash_table.h>
#include <rsys/dynamic_array.h>
+#include <rsys/dynamic_array_uint.h>
#include <rsys/dynamic_array_uchar.h>
#include <rsys/clock_time.h>
@@ -61,15 +62,23 @@ const struct cc_descriptor CC_DESCRIPTOR_NULL = CC_DESCRIPTOR_NULL__;
#define HTABLE_DATA char
#include <rsys/hash_table.h>
+/* A threshold for dot products:
+ * If ray.normal < threshold we suspect accuracy could be a problem */
+#define DOT_THRESHOLD 0.0001f
+
struct filter_ctx1 {
struct senc3d_scene* scn;
struct s3d_scene_view* view;
component_id_t origin_component;
struct darray_triangle_comp* triangles_comp;
struct darray_ptr_component_descriptor* components;
- /* Result of hit */
+ /* Tmp data used across filter calls */
double current_6volume;
+ float s;
+ /* Result of hit */
component_id_t hit_component;
+ float hit_dir[3], hit_dist;
+ struct s3d_primitive hit_prim;
};
struct filter_ctx2 {
@@ -78,9 +87,16 @@ struct filter_ctx2 {
component_id_t component;
};
+struct filter_ctx3 {
+ /* Result of hit */
+ struct darray_uint* collected_trg_ids;
+};
+
enum fctx_type {
FCTX1,
- FCTX2
+ FCTX2,
+ FCTX3,
+ FCTX4
};
struct filter_ctx {
@@ -88,6 +104,7 @@ struct filter_ctx {
union {
struct filter_ctx1 ctx1;
struct filter_ctx2 ctx2;
+ struct filter_ctx3 ctx3;
} c;
};
@@ -162,212 +179,240 @@ self_hit_filter
void* filter_data)
{
struct filter_ctx* fctx_ = ray_data;
- struct filter_ctx1* fctx;
- const struct triangle_comp* trg_comp;
- const component_id_t* hit_comp;
- float s = 0;
- enum senc3d_side hit_side;
- int i;
- double org_z, mz = -INF;
- const struct triangle_in* triangles;
- const struct triangle_in* trg = NULL;
- const union double3* vertices;
- struct cc_descriptor* const* comp_descriptors;
-
- (void)ray_dir; (void)ray_range; (void)filter_data;
- ASSERT(hit && fctx_);
-
- if(fctx_->type == FCTX2) {
- /* The filter is used to count the hits on some component along
- * an infinite ray */
- struct filter_ctx2* ctx2 = &fctx_->c.ctx2;
- ASSERT(hit->prim.prim_id
- < darray_triangle_comp_size_get(ctx2->triangles_comp));
- trg_comp = darray_triangle_comp_cdata_get(ctx2->triangles_comp);
- hit_comp = trg_comp[hit->prim.prim_id].component;
- if(hit_comp[SENC3D_FRONT] == ctx2->component
- || hit_comp[SENC3D_BACK] == ctx2->component)
- {
- ctx2->cpt++;
- }
- return 1; /* Reject to continue counting */
- }
- /* The filter is called from a point query on successive hits found from
- * ray_org, that belongs to origin_component. It can keep or reject the hit.
- * Hits are only submitted inside a certain radius from ray_org, that is
- * decreased to the hit distance for every hit that is kept.
- * At the end, the last kept hit (= the closest), determines a component to
- * which origin_component is linked. At a later stage the algorithm process
- * linked components to determine their relative inclusions.
- *
- * For each hit, the filter computes if the hit is on a component above
- * origin_component (that is with >= Z).
- * If the hit is distant (dist>0), we just keep the hit as a valid candidate,
- * but things get more tricky when dist==0 (ray_org is a vertex where some
- * other components are in contact with origin_component).
- * In this case, one of the other components can include the origin_component
- * (greater volume needed), or they can be disjoint, with (at least) ray_org
- * as a common vertex (they can also partially intersect, but this is invalid
- * and remains undetected by star enclosures). */
- ASSERT(fctx_->type == FCTX1);
- fctx = &fctx_->c.ctx1;
- comp_descriptors = darray_ptr_component_descriptor_cdata_get(fctx->components);
- trg_comp = darray_triangle_comp_cdata_get(fctx->triangles_comp);
- hit_comp = trg_comp[hit->prim.prim_id].component;
- triangles = darray_triangle_in_cdata_get(&fctx->scn->triangles_in);
- vertices = darray_position_cdata_get(&fctx->scn->vertices);
- ASSERT(hit->prim.prim_id
- < darray_triangle_comp_size_get(fctx->triangles_comp));
- ASSERT(hit->uv[0] == CLAMP(hit->uv[0], 0, 1));
- ASSERT(hit->uv[1] == CLAMP(hit->uv[1], 0, 1));
-
- /* No self hit */
- if(hit_comp[SENC3D_FRONT] == fctx->origin_component
- || hit_comp[SENC3D_BACK] == fctx->origin_component)
- return 1; /* Reject */
-
- if(hit->distance == 0) {
- /* origin_component is in contact with some other components
- * We will need further exploration to know if they should be considered */
- int n;
-
- /* If same component, process only once */
- FOR_EACH(n, 0, (hit_comp[SENC3D_FRONT] == hit_comp[SENC3D_BACK] ? 1 : 2)) {
- const enum senc3d_side sides[2] = { SENC3D_FRONT, SENC3D_BACK };
- component_id_t c = hit_comp[sides[n]];
- side_id_t side;
- double pt[3] = { 0, 0, 0 };
- float org[3], dir[3] = { 0, 0, 1 }, rg[2] = { 0, FLT_MAX };
- struct s3d_hit hit2 = S3D_HIT_NULL;
- struct filter_ctx fctx2;
- ASSERT(c < darray_ptr_component_descriptor_size_get(fctx->components));
- if(comp_descriptors[c]->is_outer_border) {
- if(fabs(comp_descriptors[c]->_6volume)
- <= fabs(comp_descriptors[fctx->origin_component]->_6volume))
- /* Component is not large enough to include origin_component */
- continue;
- } else {
- vrtx_id_t c_z_id = comp_descriptors[c]->max_z_vrtx_id;
- vrtx_id_t o_z_id = comp_descriptors[fctx->origin_component]->max_z_vrtx_id;
- ASSERT(c_z_id < darray_position_size_get(&fctx->scn->vertices));
- ASSERT(o_z_id < darray_position_size_get(&fctx->scn->vertices));
- if(vertices[c_z_id].pos.z <= vertices[o_z_id].pos.z)
- /* Component is not above origin_component */
- continue;
- }
- /* Check if component c includes origin_component; as components cannot
- * overlap, testing a single point is OK, as long as the point is not on
- * c boundary (can be on origin_component boundary, though).
- * As this case is supposed to be rare, we go for a basic algorithm */
- for(side = comp_descriptors[fctx->origin_component]->side_range.first;
- side <= comp_descriptors[fctx->origin_component]->side_range.last;
- side++)
+ (void)ray_org; (void)ray_range; (void)filter_data;
+ ASSERT(fctx_);
+
+ switch (fctx_->type) {
+ default: FATAL("Invalid");
+
+ case FCTX2: {
+ /* The filter is used to count the hits on some component along an
+ * infinite ray */
+ struct filter_ctx2* ctx2 = &fctx_->c.ctx2;
+ const struct triangle_comp* trg_comp;
+ const component_id_t* hit_comp;
+ ASSERT(hit->prim.prim_id
+ < darray_triangle_comp_size_get(ctx2->triangles_comp));
+ trg_comp = darray_triangle_comp_cdata_get(ctx2->triangles_comp);
+ hit_comp = trg_comp[hit->prim.prim_id].component;
+ if(hit_comp[SENC3D_FRONT] == ctx2->component
+ || hit_comp[SENC3D_BACK] == ctx2->component)
{
- /* Find a triangle on origin_component boundary that is not on c
- * boundary (the 2 components cannot share all their triangles) */
- trg_id_t t = TRGSIDE_2_TRG(side);
- const component_id_t* candidate_comp = trg_comp[t].component;
- if(candidate_comp[SENC3D_FRONT] != fctx->origin_component
- && candidate_comp[SENC3D_BACK] != fctx->origin_component)
- continue;
- if(candidate_comp[SENC3D_FRONT] == c
- || candidate_comp[SENC3D_BACK] == c)
- continue;
- /* This triangle is OK */
- trg = triangles + t;
- break;
- }
- ASSERT(trg != NULL);
- /* Any point on trg not on an edge can do the trick: use the barycenter */
- FOR_EACH(i, 0, 3) {
- vrtx_id_t v = trg->vertice_id[i];
- ASSERT(v < darray_position_size_get(&fctx->scn->vertices));
- d3_add(pt, pt, vertices[v].vec);
- }
- d3_divd(pt, pt, 3);
- f3_set_d3(org, pt);
- /* Trace a ray and count intersections with components of trg */
- fctx2.type = FCTX2;
- fctx2.c.ctx2.triangles_comp = fctx->triangles_comp;
- fctx2.c.ctx2.cpt = 0;
- fctx2.c.ctx2.component = c;
- S3D(scene_view_trace_ray(fctx->view, org, dir, rg, &fctx2, &hit2));
- ASSERT(S3D_HIT_NONE(&hit2)); /* The ray is supposed to go to infinity */
- /* origin_component is linked_to an outer component if cpt is odd,
- * linked_to an inner component if cpt is even */
- if(comp_descriptors[c]->is_outer_border == (fctx2.c.ctx2.cpt % 2)) {
- double v = fabs(comp_descriptors[c]->_6volume);
- /* If origin_component is inside several components, the one we are
- * looking for is the smallest one */
- if(v >= fctx->current_6volume) continue;
- fctx->hit_component = fctx2.c.ctx2.component;
- fctx->current_6volume = v;
- /* Continue searching for a smaller component that includes
- * origin_component */
+ ctx2->cpt++;
}
+ return 1; /* Reject to continue counting */
}
- return 1; /* Reject */
- }
-
- /* Reject hits with < Z */
- ASSERT(hit->prim.prim_id <
- darray_triangle_in_size_get(&fctx->scn->triangles_in));
- trg = triangles + hit->prim.prim_id;
- /* Check if the hit is above ray_org (hit[2] > ray_org[2])
- * As we cannot rely on numerical accuracy when computing hit positions,
- * we use the triangle vertices to check if some part of the hit triangle
- * is above ray_org */
- FOR_EACH(i, 0, 3) {
- vrtx_id_t v = trg->vertice_id[i];
- const union double3* p = vertices + v;
- ASSERT(v < darray_position_size_get(&fctx->scn->vertices));
- if(i == 0 || mz < p->pos.z) mz = p->pos.z;
- }
- /* Don't use org[2] as, being float, it would lead to a float VS double
- * comparison that causes accuracy problems. */
- org_z = vertices[comp_descriptors[fctx->origin_component]->max_z_vrtx_id].pos.z;
- if(mz <= org_z)
- return 1; /* Hit triangle is below ray_org: reject */
-
- if(hit_comp[SENC3D_FRONT] == hit_comp[SENC3D_BACK]) {
- /* Easy case and hit component is known */
- fctx->hit_component = hit_comp[SENC3D_FRONT];
- return 0; /* Keep */
- }
- /* Compute hit side by using the vertex with best accuracy */
- FOR_EACH(i, 0, 3) {
- float tmps, tmp[3], dir[3];
- vrtx_id_t v = trg->vertice_id[i];
- const union double3* p = vertices + v;
- ASSERT(v < darray_position_size_get(&fctx->scn->vertices));
- f3_sub(dir, f3_set_d3(tmp, p->vec), ray_org);
- tmps = f3_dot(dir, hit->normal);
- if(i == 0 || fabsf(s) < fabsf(tmps)) s = tmps;
- }
+ case FCTX3: {
+ /* The filter is used to posprocess FCTX1 cases ending with a nearly
+ * parallel triangle that does not allow to select reliably a side and its
+ * attached component information using a closest point request.
+ *
+ * This filter will collect any triangle in the search area (around the
+ * point returned at step FCTX1). These triangles will be further
+ * processed by FCTX4 filter with a raytracing request. */
+ struct filter_ctx3* ctx = &fctx_->c.ctx3;
+ res_T res = darray_uint_push_back(ctx->collected_trg_ids, &hit->prim.prim_id);
+ CHK(res == RES_OK);
+ return 1; /* Reject do avoid decreasing the search radius */
+ }
- /* We cannot know which side to consider if s==0.
- * As hit distance > 0 and the 2 sides belong to 2 different components,
- * another candidate must selected afterwards (can be at greater distance,
- * disallowing to restrict the search distance here) */
- if(s == 0) {
- fctx->hit_component = COMPONENT_NULL__;
- return 1; /* Reject */
- }
+ case FCTX4:
+ /* The filter is used to posprocess FCTX3 triangles using a single ray.
+ *
+ * No specific processing is due here. */
+ return 0; /* Keep */
+
+ case FCTX1: {
+ /* This filter is called from a closest point query from a point belonging
+ * to origin_component. The returned hit is used to determine a component
+ * to which origin_component is linked. At a later stage the algorithm
+ * process linked components to determine their relative inclusions.
+ *
+ * For each hit, the filter computes if the hit is on a component above
+ * origin_component (that is with >= Z).
+ * If the hit is distant (dist>0), we just keep the hit as a valid candidate,
+ * but things get more tricky when dist==0 (ray_org is a vertex where some
+ * other components can be in contact with origin_component).
+ * In this case, one of the other components can include the origin_component
+ * (greater volume needed), or they can be disjoint, with (at least) ray_org
+ * as a common vertex (they can also partially intersect, but this is invalid
+ * and remains undetected by star enclosures). */
+ struct filter_ctx1* ctx = &fctx_->c.ctx1;
+ struct cc_descriptor* const*
+ comp_descriptors = darray_ptr_component_descriptor_cdata_get(ctx->components);
+ const struct triangle_comp*
+ trg_comp = darray_triangle_comp_cdata_get(ctx->triangles_comp);
+ const component_id_t*
+ hit_comp = trg_comp[hit->prim.prim_id].component;
+ const struct triangle_in*
+ triangles = darray_triangle_in_cdata_get(&ctx->scn->triangles_in);
+ const union double3* vertices = darray_position_cdata_get(&ctx->scn->vertices);
+ const struct triangle_in* trg = NULL;
+ double org_z, mz = -INF;
+ enum senc3d_side hit_side;
+ float s = 0, hit_normal[3], rdir[3];
+ int i;
- /* Determine which side was hit */
- hit_side =
- ((s < 0) /* Facing geometrical normal of hit */
- == ((fctx->scn->convention & SENC3D_CONVENTION_NORMAL_FRONT) != 0))
- /* Warning: following Embree 2 convention for geometrical normals,
- * the Star3D hit normal is left-handed while star-enclosures-3d uses
- * right-handed convention */
- ? SENC3D_BACK : SENC3D_FRONT;
+ ASSERT(hit->prim.prim_id
+ < darray_triangle_comp_size_get(ctx->triangles_comp));
+ ASSERT(hit->uv[0] == CLAMP(hit->uv[0], 0, 1));
+ ASSERT(hit->uv[1] == CLAMP(hit->uv[1], 0, 1));
+
+ /* No self hit */
+ if(hit_comp[SENC3D_FRONT] == ctx->origin_component
+ || hit_comp[SENC3D_BACK] == ctx->origin_component)
+ return 1; /* Reject */
+
+ if(hit->distance == 0) {
+ /* origin_component is in contact with some other components
+ * We will need further exploration to know if they should be considered */
+ int n;
+
+ /* If same component, process only once */
+ FOR_EACH(n, 0, (hit_comp[SENC3D_FRONT] == hit_comp[SENC3D_BACK] ? 1 : 2)) {
+ const enum senc3d_side sides[2] = { SENC3D_FRONT, SENC3D_BACK };
+ component_id_t c = hit_comp[sides[n]];
+ side_id_t side;
+ double pt[3] = { 0, 0, 0 };
+ float org[3], dir[3] = { 0, 0, 1 }, rg[2] = { 0, FLT_MAX };
+ struct s3d_hit hit2 = S3D_HIT_NULL;
+ struct filter_ctx fctx2;
+ ASSERT(c < darray_ptr_component_descriptor_size_get(ctx->components));
+ if(comp_descriptors[c]->is_outer_border) {
+ if(fabs(comp_descriptors[c]->_6volume)
+ <= fabs(comp_descriptors[ctx->origin_component]->_6volume))
+ /* Component is not large enough to include origin_component */
+ continue;
+ } else {
+ vrtx_id_t c_z_id = comp_descriptors[c]->max_z_vrtx_id;
+ vrtx_id_t o_z_id = comp_descriptors[ctx->origin_component]->max_z_vrtx_id;
+ ASSERT(c_z_id < darray_position_size_get(&ctx->scn->vertices));
+ ASSERT(o_z_id < darray_position_size_get(&ctx->scn->vertices));
+ if(vertices[c_z_id].pos.z <= vertices[o_z_id].pos.z)
+ /* Component is not above origin_component */
+ continue;
+ }
+ /* Check if component c includes origin_component; as components cannot
+ * overlap, testing a single point is OK, as long as the point is not on
+ * c boundary (can be on origin_component boundary, though).
+ * As this case is supposed to be rare, we go for a basic algorithm */
+ for(side = comp_descriptors[ctx->origin_component]->side_range.first;
+ side <= comp_descriptors[ctx->origin_component]->side_range.last;
+ side++)
+ {
+ /* Find a triangle on origin_component boundary that is not on c
+ * boundary (the 2 components cannot share all their triangles) */
+ trg_id_t t = TRGSIDE_2_TRG(side);
+ const component_id_t* candidate_comp = trg_comp[t].component;
+ if(candidate_comp[SENC3D_FRONT] != ctx->origin_component
+ && candidate_comp[SENC3D_BACK] != ctx->origin_component)
+ continue;
+ if(candidate_comp[SENC3D_FRONT] == c
+ || candidate_comp[SENC3D_BACK] == c)
+ continue;
+ /* This triangle is OK */
+ trg = triangles + t;
+ break;
+ }
+ ASSERT(trg != NULL);
+ /* Any point on trg not on an edge can do the trick: use the barycenter */
+ FOR_EACH(i, 0, 3) {
+ vrtx_id_t v = trg->vertice_id[i];
+ ASSERT(v < darray_position_size_get(&ctx->scn->vertices));
+ d3_add(pt, pt, vertices[v].vec);
+ }
+ d3_divd(pt, pt, 3);
+ f3_set_d3(org, pt);
+ /* Trace a ray and count intersections with components of trg */
+ fctx2.type = FCTX2;
+ fctx2.c.ctx2.triangles_comp = ctx->triangles_comp;
+ fctx2.c.ctx2.cpt = 0;
+ fctx2.c.ctx2.component = c;
+ S3D(scene_view_trace_ray(ctx->view, org, dir, rg, &fctx2, &hit2));
+ ASSERT(S3D_HIT_NONE(&hit2)); /* The ray is supposed to go to infinity */
+ /* origin_component is linked_to an outer component if cpt is odd,
+ * linked_to an inner component if cpt is even */
+ if(comp_descriptors[c]->is_outer_border == (fctx2.c.ctx2.cpt % 2)) {
+ double v = fabs(comp_descriptors[c]->_6volume);
+ /* If origin_component is inside several components, the one we are
+ * looking for is the smallest one */
+ if(v >= ctx->current_6volume) continue;
+ ctx->hit_component = fctx2.c.ctx2.component;
+ ctx->current_6volume = v;
+ /* Continue searching for a smaller component that includes
+ * origin_component */
+ }
+ }
+ return 1; /* Reject */
+ }
- fctx->hit_component = hit_comp[hit_side];
+ /* Reject hits with < Z */
+ ASSERT(hit->prim.prim_id <
+ darray_triangle_in_size_get(&ctx->scn->triangles_in));
+ trg = triangles + hit->prim.prim_id;
+ /* Check if the hit is above ray_org (hit[2] > ray_org[2])
+ * As we cannot rely on numerical accuracy when computing hit positions,
+ * we use the triangle vertices to check if some part of the hit triangle
+ * is above ray_org */
+ FOR_EACH(i, 0, 3) {
+ vrtx_id_t v = trg->vertice_id[i];
+ const union double3* p = vertices + v;
+ ASSERT(v < darray_position_size_get(&ctx->scn->vertices));
+ if(i == 0 || mz < p->pos.z) mz = p->pos.z;
+ }
+ /* Don't use org[2] as, being float, it would lead to a float VS double
+ * comparison that causes accuracy problems. */
+ org_z = vertices[comp_descriptors[ctx->origin_component]->max_z_vrtx_id].pos.z;
+ if(mz <= org_z)
+ return 1; /* Hit triangle is below ray_org: reject */
+
+ if(hit_comp[SENC3D_FRONT] == hit_comp[SENC3D_BACK]) {
+ /* Easy case and hit component is known */
+ ctx->hit_component = hit_comp[SENC3D_FRONT];
+ ctx->s = 1;
+ f3_set(ctx->hit_dir, ray_dir);
+ ctx->hit_dist = hit->distance;
+ ctx->hit_prim = hit->prim;
+ return 0; /* Keep */
+ }
- return 0; /* Keep */
+ /* Compute hit side */
+ /* For s to be comparable, vectors must be normailzed */
+ f3_normalize(hit_normal, hit->normal);
+ f3_normalize(rdir, ray_dir);
+ s = f3_dot(rdir, hit_normal);
+
+ if(fabsf(s) < DOT_THRESHOLD) {
+ /* We cannot know for sure which side to consider */
+ if(ctx->hit_dist == hit->distance && ctx->s > s) {
+ /* Same distance with worse s: keep the previous hit */
+ return 1; /* Reject */
+ }
+ ctx->hit_component = COMPONENT_NULL__;
+ ctx->s = s;
+ f3_set(ctx->hit_dir, ray_dir);
+ ctx->hit_dist = hit->distance;
+ ctx->hit_prim = hit->prim;
+ return 0; /* Keep */
+ }
+ /* Determine which side was hit */
+ hit_side =
+ ((s < 0) /* Facing geometrical normal of hit */
+ == ((ctx->scn->convention & SENC3D_CONVENTION_NORMAL_FRONT) != 0))
+ /* Warning: following Embree 2 convention for geometrical normals,
+ * the Star3D hit normal is left-handed while star-enclosures-3d uses
+ * right-handed convention */
+ ? SENC3D_BACK : SENC3D_FRONT;
+ ctx->hit_component = hit_comp[hit_side];
+ ctx->s = s;
+ f3_set(ctx->hit_dir, ray_dir);
+ ctx->hit_dist = hit->distance;
+ ctx->hit_prim = hit->prim;
+ return 0; /* Keep */
+ }
+ }
}
static void
@@ -425,7 +470,7 @@ extract_connex_components
const struct side_range* media_use
= darray_side_range_cdata_get(&scn->media_use);
FOR_EACH(s, 0, 2) {
- const side_id_t side = TRGIDxSIDE_2_TRGSIDE(t_, s);
+ const side_id_t side = TRGIDxSIDE_2_TRGSIDE(t_, (enum senc3d_side)s);
medium_id_t medium = trg_in->medium[s];
m_idx = medium_id_2_medium_idx(medium);
ASSERT(media_use[m_idx].first <= side && side
@@ -893,8 +938,10 @@ group_connex_components
size_t tmp;
component_id_t cc_count;
int64_t ccc;
- struct filter_ctx fctx;
+ struct filter_ctx ctx1, ctx3, ctx4;
float lower[3], upper[3];
+ struct darray_uint collected_triangles;
+ int collected_initialized = 0;
ASSERT(scn && triangles_comp && connex_components
&& s3d_view && next_enclosure_id && res);
@@ -914,14 +961,13 @@ group_connex_components
ASSERT(tmp <= COMPONENT_MAX__);
cc_count = (component_id_t)tmp;
positions = darray_position_cdata_get(&scn->vertices);
- fctx.type = FCTX1;
- fctx.c.ctx1.scn = scn;
- fctx.c.ctx1.view = s3d_view;
- fctx.c.ctx1.triangles_comp = triangles_comp;
- fctx.c.ctx1.components = connex_components;
+ ctx1.type = FCTX1;
+ ctx1.c.ctx1.scn = scn;
+ ctx1.c.ctx1.view = s3d_view;
+ ctx1.c.ctx1.triangles_comp = triangles_comp;
+ ctx1.c.ctx1.components = connex_components;
*res = s3d_scene_view_get_aabb(s3d_view, lower, upper);
if(*res != RES_OK) goto end;
- /* Cast rays to find links between connex components */
#pragma omp for schedule(dynamic)
for(ccc = 0; ccc < (int64_t)cc_count; ccc++) {
res_T tmp_res = RES_OK;
@@ -950,9 +996,10 @@ group_connex_components
f3_set_d3(origin, max_vrtx);
/* Self-hit data: self hit if hit this component "on the other side" */
- fctx.c.ctx1.origin_component = cc->cc_id;
- fctx.c.ctx1.current_6volume = DBL_MAX;
- fctx.c.ctx1.hit_component = COMPONENT_NULL__;
+ ctx1.c.ctx1.origin_component = cc->cc_id;
+ ctx1.c.ctx1.current_6volume = DBL_MAX;
+ ctx1.c.ctx1.hit_dist = FLT_MAX;
+ ctx1.c.ctx1.hit_component = COMPONENT_NULL__;
/* Limit search radius. Only upwards moves (+Z) will be considered.
* Use a radius allowing to reach the closest top vertex of scene's AABB */
FOR_EACH(i, 0, 2) {
@@ -962,18 +1009,131 @@ group_connex_components
ASSERT(lower[2] <= origin[2] && origin[2] <= upper[2]);
rrr[2] = upper[2] - origin[2];
r = f3_len(rrr) + FLT_EPSILON; /* Ensure r > 0 */
- tmp_res = s3d_scene_view_closest_point(s3d_view, origin, r, &fctx, &hit);
+ /* Cast a sphere to find links between connex components */
+ tmp_res = s3d_scene_view_closest_point(s3d_view, origin, r, &ctx1, &hit);
if(tmp_res != RES_OK) {
*res = tmp_res;
continue;
}
/* If no hit, the component is facing an infinite medium */
- if(fctx.c.ctx1.hit_component == COMPONENT_NULL__) {
+ if(ctx1.c.ctx1.hit_dist == FLT_MAX) {
cc->cc_group_root = CC_GROUP_ROOT_INFINITE;
cc->enclosure_id = 0;
+ }
+ else if(ctx1.c.ctx1.hit_component == COMPONENT_NULL__) {
+ /* The selected triangle was nearly parallel to the line of sight:
+ * FRONT/BACK discrimination was not reliable enough and should be done
+ * differently. */
+ const struct triangle_in*
+ triangles = darray_triangle_in_cdata_get(&scn->triangles_in);
+ const union double3* vertices = darray_position_cdata_get(&scn->vertices);
+ struct s3d_hit tmp_hit = S3D_HIT_NULL;
+ const struct triangle_comp*
+ trg_comp = darray_triangle_comp_cdata_get(triangles_comp);
+ component_id_t hit_comp;
+ float hit_pt[3];
+ float short_range;
+ const unsigned* ids;
+ size_t n;
+ int found_conclusive = 0;
+ if(!collected_initialized) {
+ darray_uint_init(scn->dev->allocator, &collected_triangles);
+ collected_initialized = 1;
+ }
+ /* First step is to collect triangles around the hit previously returned */
+ ctx3.type = FCTX3;
+ darray_uint_clear(&collected_triangles);
+ ctx3.c.ctx3.collected_trg_ids = &collected_triangles;
+ /* Check that |ray_dir| == ray_dist i.e. no need to use ray_dir*ray_dist */
+ ASSERT(f3_len(ctx1.c.ctx1.hit_dir) == ctx1.c.ctx1.hit_dist);
+ f3_add(hit_pt, origin, ctx1.c.ctx1.hit_dir);
+ short_range = MMAX(1e-3f, 1e-3f * f3_len(hit_pt));
+ tmp_res = s3d_scene_view_closest_point(s3d_view, hit_pt, short_range,
+ &ctx3, &tmp_hit);
+ if(tmp_res != RES_OK || darray_uint_size_get(&collected_triangles) == 0) {
+ log_err(scn->dev, LIB_NAME": %s: %d: collect failed.\n", FUNC_NAME,
+ __LINE__ );
+ *res = (tmp_res != RES_OK) ? tmp_res : RES_BAD_OP;
+ continue;
+ }
+ /* Second step is to trace rays towards these triangles to find the
+ * closest geometry from origin whose normal allows to discrimine reliably
+ * which side is seen. Stop after a first conclusive triangle. */
+ ids = darray_uint_cdata_get(&collected_triangles);
+ for(n = 0; n < darray_uint_size_get(&collected_triangles); n++) {
+ const struct triangle_in* trg = triangles + ids[n];
+ float rdir[3], range[2], bary[3];
+ double pt[3] = { 0, 0, 0 };
+ int fst_hit_found = 0;
+ ctx4.type = FCTX4;
+ /* Any point on trg can do the trick: use the barycenter */
+ FOR_EACH(i, 0, 3) {
+ vrtx_id_t v = trg->vertice_id[i];
+ ASSERT(v < darray_position_size_get(&scn->vertices));
+ d3_add(pt, pt, vertices[v].vec);
+ }
+ f3_set_d3(bary, pt);
+ f3_divf(bary, bary, 3);
+ f3_sub(rdir, bary, origin);
+ f3_normalize(rdir, rdir);
+ /* Search starts at half the distance of the closest geometry to avoid
+ * self hits */
+ range[0] = ctx1.c.ctx1.hit_dist * 0.5f;
+ range[1] = FLT_MAX;
+ while(1) {
+ S3D(scene_view_trace_ray(s3d_view, origin, rdir, range, &ctx4, &tmp_hit));
+ if(S3D_HIT_NONE(&tmp_hit)) {
+ /* Unexpected => next candidate */
+ break;
+ }
+ if(!fst_hit_found) {
+ /* Continue only if conclusive */
+ float normal[3], s;
+ enum senc3d_side hit_side;
+ f3_normalize(normal, tmp_hit.normal);
+ s = f3_dot(rdir, normal);
+ if(fabsf(s) < DOT_THRESHOLD) {
+ /* Unconclusive */
+ break;
+ }
+ /* Determine which side was hit */
+ hit_side =
+ ((s < 0) /* Facing geometrical normal of hit */
+ == ((scn->convention & SENC3D_CONVENTION_NORMAL_FRONT) != 0))
+ /* Warning: following Embree 2 convention for geometrical normals,
+ * the Star3D hit normal is left-handed while star-enclosures-3d uses
+ * right-handed convention */
+ ? SENC3D_BACK : SENC3D_FRONT;
+ /* Keep it */
+ hit_comp = trg_comp[tmp_hit.prim.prim_id].component[hit_side];
+ fst_hit_found = 1;
+ }
+ if(!S3D_PRIMITIVE_EQ(&ctx1.c.ctx1.hit_prim, &tmp_hit.prim)) {
+ /* Not a problem as long as the ray eventualy hits the barycenter */
+ range[0] = nextafterf(tmp_hit.distance, FLT_MAX);
+ continue; /* raytrace again to visit the end of the segment */
+ }
+ /* At the barycenter */
+ break;
+ }
+ if(fst_hit_found) {
+ found_conclusive = 1;
+ break; /* Don't need more */
+ }
+ }
+ if(!found_conclusive) {
+ /* Could try something else; now just report a failure */
+ log_err(scn->dev, LIB_NAME": %s: %d: search failed.\n", FUNC_NAME,
+ __LINE__ );
+ *res = RES_BAD_OP;
+ continue;
+ }
+ /* We finally got it! Group this component */
+ cc->cc_group_root = hit_comp;
+ ASSERT(cc->cc_group_root < cc_count);
} else {
/* If hit, group this component */
- cc->cc_group_root = fctx.c.ctx1.hit_component;
+ cc->cc_group_root = ctx1.c.ctx1.hit_component;
ASSERT(cc->cc_group_root < cc_count);
}
}
@@ -1031,6 +1191,7 @@ group_connex_components
}
/* Implicit barrier here */
end:
+ if(collected_initialized) darray_uint_release(&collected_triangles);
return;
}
@@ -1179,7 +1340,7 @@ collect_and_link_neighbours
/* Compute rotation angle around common edge */
d3_sub(edge, vertices[v2].vec, vertices[v0].vec);
d33_muld3(edge, basis, edge);
- ASSERT(d3_len(edge) && (edge[0] || edge[1]));
+ ASSERT(d3_len(edge) != 0 && (edge[0] != 0 || edge[1] != 0));
neighbour_info->angle = atan2(edge[1], edge[0]); /* in ]-pi + pi]*/
if(is_reversed)
d3(n.vec, +edge[1], -edge[0], 0);
