star-enclosures-3d

Extract enclosures from 3D geometry
git clone git://git.meso-star.fr/star-enclosures-3d.git
Log | Files | Refs | README | LICENSE

commit 2d32ffff33b6d235f1db514771cdea13a19f07d1
parent 906ed56bcfc2e8c105602b7b73dd6d9b3f68c57f
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Thu,  7 May 2020 09:55:55 +0200

Add area and volume to enclosure headers

Diffstat:
Msrc/senc3d.h | 14++++++++++++--
Msrc/senc3d_enclosure_data.h | 4+++-
Msrc/senc3d_scene_analyze.c | 139+++++++++++++++++++++++++++++++++++++++++++------------------------------------
Msrc/senc3d_scene_analyze_c.h | 6+++++-
Msrc/test_senc3d_enclosure.c | 4++++
Msrc/test_senc3d_inconsistant_cube.c | 20++++++++++++++++----
6 files changed, 115 insertions(+), 72 deletions(-)

diff --git a/src/senc3d.h b/src/senc3d.h @@ -76,9 +76,9 @@ enum senc3d_side { struct senc3d_enclosure_header { /* The ID of the enclosure; 0, 1, ... */ unsigned enclosure_id; - /* Number of triangles; a triangle can be accounted for twice, once by side*/ + /* Number of triangles; a triangle can be counted twice, once by side*/ unsigned primitives_count; - /* Number of unique triangles; a triangle cannot be accounted for twice */ + /* Number of unique triangles; a triangle cannot be counted twice */ unsigned unique_primitives_count; /* Number of vertices */ unsigned vertices_count; @@ -88,6 +88,16 @@ struct senc3d_enclosure_header { /* Is the enclosure open/infinite? * Only the outermost enclosure is infinite. */ int is_infinite; + /* The volume of the enclosure, in m^3. + * For the outermost enclosure, that goes to infinity, the volume is negative + * and is the volume substracted from an hypothetical surrounding object. + * If the two sides of a triangle are part of the enclosure, the triangle is + * counted as 0. */ + double volume; + /* The area of the enclosure, in m^2. + * If the two sides of a triangle are part of the enclosure, the triangle is + * counted twice. */ + double area; }; /* We consider the geometrical normal Ng to a triangle V0 V1 V2 diff --git a/src/senc3d_enclosure_data.h b/src/senc3d_enclosure_data.h @@ -66,7 +66,9 @@ init_header(struct senc3d_enclosure_header* header) header->unique_primitives_count = 0; header->vertices_count = 0; header->enclosed_media_count = 0; - header->is_infinite = CHAR_MAX; + header->is_infinite = INT_MAX; + header->volume = 0; + header->area = 0; } #define DARRAY_NAME media diff --git a/src/senc3d_scene_analyze.c b/src/senc3d_scene_analyze.c @@ -35,7 +35,7 @@ #include <stdlib.h> #define CC_DESCRIPTOR_NULL__ {\ - CHAR_MAX, VRTX_NULL__, 0,\ + 0, 0, INT_MAX, VRTX_NULL__, 0,\ CC_ID_NONE, CC_GROUP_ROOT_NONE, ENCLOSURE_NULL__,\ { TRG_NULL__, 0},\ NULL\ @@ -188,7 +188,6 @@ extract_connex_components struct mem_allocator* alloc; int64_t m_idx; /* OpenMP requires a signed type for the for loop variable */ struct darray_side_id stack; - struct darray_side_id ids_of_sides_around_max_z_vertex; const union double3* positions; const struct triangle_tmp* triangles_tmp; struct triangle_comp* triangles_comp; @@ -207,7 +206,6 @@ extract_connex_components triangles_tmp = darray_triangle_tmp_cdata_get(triangles_tmp_array); triangles_comp = darray_triangle_comp_data_get(triangles_comp_array); darray_side_id_init(alloc, &stack); - darray_side_id_init(alloc, &ids_of_sides_around_max_z_vertex); darray_side_id_init(alloc, &current_component); processed = MEM_CALLOC(alloc, scn->ntris, sizeof(uchar)); if(!processed) { @@ -304,25 +302,22 @@ extract_connex_components enum side_flag crt_side_flag = TRGSIDE_2_SIDEFLAG(crt_side_id); struct trgside* crt_side = trgsides + crt_side_id; const trg_id_t crt_trg_id = TRGSIDE_2_TRG(crt_side_id); - const struct triangle_in* trg_in = - darray_triangle_in_cdata_get(&scn->triangles_in) + crt_trg_id; uchar* trg_used = processed + crt_trg_id; const struct triangle_tmp* const trg_tmp = triangles_tmp + crt_trg_id; ASSERT(crt_trg_id < scn->ntris); if(*p_res != RES_OK) break; - /* Record Zmax information - * Keep track of the appropriate vertex of the component in order - * to cast a ray at the component grouping step of the algorithm. - * The most appropriate vertex is (the) one with the greater Z - * coordinate. */ + /* Record max_z information + * Keep track of the appropriate vertex of the component in order to + * start upwards closest point query at the component grouping step of + * the algorithm. The most appropriate vertex is (the) one with the + * greater Z coordinate. */ if(max_z < trg_tmp->max_z) { /* New best vertex */ + const struct triangle_in* trg_in = + darray_triangle_in_cdata_get(&scn->triangles_in) + crt_trg_id; max_z = trg_tmp->max_z; - /* New vertex: reset list of sides */ - darray_side_id_clear(&ids_of_sides_around_max_z_vertex); - /* Select a vertex with z == max_z */ FOR_EACH(i, 0, 3) { if(max_z == positions[trg_in->vertice_id[i]].pos.z) { @@ -331,19 +326,6 @@ extract_connex_components } } ASSERT(i < 3); /* Found one */ - /* List of sides using the vertex */ - OK2(darray_side_id_push_back(&ids_of_sides_around_max_z_vertex, - &crt_side_id)); - } else if(max_z == trg_tmp->max_z) { - /* Does this triangle use the currently selected max_z vertex? */ - FOR_EACH(i, 0, 3) { - if(max_z_vrtx_id == trg_in->vertice_id[i]) { - /* List of sides using the vertex */ - OK2(darray_side_id_push_back(&ids_of_sides_around_max_z_vertex, - &crt_side_id)); - break; - } - } } /* Record crt_side both as component and triangle level */ @@ -434,53 +416,80 @@ extract_connex_components {STATIC_ASSERT(sizeof(cc->cc_id) >= 4, Cannot_write_IDs_sync_free);} ASSERT(IS_ALIGNED(triangles_comp->component, sizeof(cc->cc_id))); FOR_EACH(ii, 0, sz) { - const side_id_t s = darray_side_id_cdata_get(&current_component)[ii]; - trg_id_t tid = TRGSIDE_2_TRG(s); - enum senc3d_side sid = TRGSIDE_2_SIDE(s); - component_id_t* cmp = triangles_comp[tid].component; - ASSERT(cmp[sid] == COMPONENT_NULL__); - ASSERT(medium_id_2_medium_idx(trgsides[s].medium) + const side_id_t s_id = darray_side_id_cdata_get(&current_component)[ii]; + trg_id_t t_id = TRGSIDE_2_TRG(s_id); + enum senc3d_side side = TRGSIDE_2_SIDE(s_id); + component_id_t* cmp = triangles_comp[t_id].component; + ASSERT(cmp[side] == COMPONENT_NULL__); + ASSERT(medium_id_2_medium_idx(trgsides[s_id].medium) >= medium_id_2_medium_idx(medium)); - cmp[sid] = cc->cc_id; + cmp[side] = cc->cc_id; } - - /* Determine if this component can be an inner part inside another - * component (creating a hole) as only these components will need - * to search for their possible outer component when grouping - * components to create encosures */ + /* Compute component area and volume, and record information on the + * max_z side of the component to help find out if the component is + * inner or outer */ max_nz = 0; max_nz_side_id = SIDE_NULL__; - sz = darray_side_id_size_get(&ids_of_sides_around_max_z_vertex); - ASSERT(sz > 0); FOR_EACH(ii, 0, sz) { - const side_id_t side_id = - darray_side_id_cdata_get(&ids_of_sides_around_max_z_vertex)[ii]; - const trg_id_t trg_id = TRGSIDE_2_TRG(side_id); - enum senc3d_side s = TRGSIDE_2_SIDE(side_id); + const side_id_t s_id = darray_side_id_cdata_get(&current_component)[ii]; + trg_id_t t_id = TRGSIDE_2_TRG(s_id); + enum senc3d_side side = TRGSIDE_2_SIDE(s_id); + struct triangle_comp* trg_comp = triangles_comp + t_id; + const struct triangle_tmp* const trg_tmp = triangles_tmp + t_id; const struct triangle_in* trg_in = - darray_triangle_in_cdata_get(&scn->triangles_in) + trg_id; - const struct triangle_comp* trg_comp = triangles_comp + trg_id; + darray_triangle_in_cdata_get(&scn->triangles_in) + t_id; const union double3* vertices = darray_position_cdata_get(&scn->vertices); double edge0[3], edge1[3], normal[3], norm; - - ASSERT(trg_comp->component[s] == cc->cc_id); (void)s; - d3_sub(edge0, vertices[trg_in->vertice_id[1]].vec, - vertices[trg_in->vertice_id[0]].vec); - d3_sub(edge1, vertices[trg_in->vertice_id[2]].vec, - vertices[trg_in->vertice_id[0]].vec); + const double* v0 = vertices[trg_in->vertice_id[0]].vec; + const double* v1 = vertices[trg_in->vertice_id[1]].vec; + const double* v2 = vertices[trg_in->vertice_id[2]].vec; + int is_2sided = (trg_comp->component[SENC3D_FRONT] + == trg_comp->component[SENC3D_BACK]); + + /* Compute component area and volume */ + d3_sub(edge0, v1, v0); + d3_sub(edge1, v2, v0); d3_cross(normal, edge0, edge1); norm = d3_normalize(normal, normal); - ASSERT(norm); (void)norm; + ASSERT(norm); + cc->_2area += norm; + + if(!is_2sided) { + /* Build a tetrahedron whose base is the triangle and whose apex is + * the coordinate system's origin. If the 2 sides of the triangle + * are part of the component, the contribution of the triangle + * must be 0. To achieve this, we just skip both sides */ + double _2base = norm; /* 2 * base area */ + double h = -d3_dot(normal, v0); /* Height of the tetrahedron */ + if((trg_in->medium[SENC3D_FRONT] == medium) + == (scn->convention & SENC3D_CONVENTION_NORMAL_FRONT)) + cc->_6volume += (h * _2base); + else + cc->_6volume -= (h * _2base); + } - if(fabs(max_nz) < fabs(normal[2])) { - max_nz_side_id = side_id; - max_nz = normal[2]; - max_z_is_2sided = (trg_comp->component[SENC3D_FRONT] - == trg_comp->component[SENC3D_BACK]); + ASSERT(trg_comp->component[side] != COMPONENT_NULL__); (void)side; + if(trg_tmp->max_z == max_z) { + int i; + /* candidate to define the max_z normal */ + FOR_EACH(i, 0, 3) { + if(cc->max_z_vrtx_id == trg_in->vertice_id[i]) { + if(fabs(normal[2]) > fabs(max_nz)) { + max_nz_side_id = s_id; + max_nz = normal[2]; + max_z_is_2sided = is_2sided; + break; + } + } + } } } - /* The inner/outer property comes from the normal orientation of the + /* Determine if this component can be an inner part inside another + * component (creating a hole) as only these components will need + * to search for their possible outer component when grouping + * components to create enclosures. + * The inner/outer property comes from the normal orientation of the * triangle on top of the component, this triangle being the one whose * |Nz| is maximal. If this triangle has its 2 sides in the component, * the component is inner */ @@ -526,7 +535,6 @@ extract_connex_components MEM_RM(alloc, current_media); darray_side_id_release(&current_component); darray_side_id_release(&stack); - darray_side_id_release(&ids_of_sides_around_max_z_vertex); /* The first thread here creates the s3d view */ #pragma omp single nowait @@ -994,7 +1002,7 @@ compare_enclosures const struct enclosure_data* e1 = ptr1; const struct enclosure_data* e2 = ptr2; ASSERT(e1->side_range.first != e2->side_range.first); - return (e1->side_range.first - e2->side_range.first); + return (int)e1->side_range.first - (int)e2->side_range.first; } static void @@ -1060,10 +1068,13 @@ build_result enclosure_id_t rank = enclosures[e].header.enclosure_id; ordered_ids[rank] = e; } - /* Rewrite cc_descriptors */ FOR_EACH(d, 0, cc_count) { - cc_descriptors[d]->enclosure_id = - ordered_ids[cc_descriptors[d]->enclosure_id]; + enclosure_id_t new_id = ordered_ids[cc_descriptors[d]->enclosure_id]; + /* Update the enclosure ID in cc_descriptor */ + cc_descriptors[d]->enclosure_id = new_id; + /* Sum up areas and volumes of components int oenclosures */ + enclosures[new_id].header.area += cc_descriptors[d]->_2area * 0.5; + enclosures[new_id].header.volume += cc_descriptors[d]->_6volume / 6; } single_err: (void)0; }/* Implicit barrier here. */ diff --git a/src/senc3d_scene_analyze_c.h b/src/senc3d_scene_analyze_c.h @@ -43,8 +43,12 @@ init_trgside(struct mem_allocator* alloc, struct trgside* data) * Also keeps the maximum z info of the component * and the list of media found */ struct cc_descriptor { + /* Twice the area of the component */ + double _2area; + /* Six times the signed volume of the component */ + double _6volume; /* Does this component is an outer border of an enclosure or an inner one? */ - char is_outer_border; + int is_outer_border; vrtx_id_t max_z_vrtx_id; /* id of the vrtx with max z value */ side_id_t side_count; /* Used when grouping components to form enclosures */ diff --git a/src/test_senc3d_enclosure.c b/src/test_senc3d_enclosure.c @@ -39,6 +39,7 @@ test(const int convention) unsigned gid; enum senc3d_side side; double vrtx[3]; + double ext_v = 0; struct context ctx = CONTEXT_NULL__; unsigned i, n, t, count; int conv; @@ -136,6 +137,8 @@ test(const int convention) CHK(header.unique_primitives_count == ntriangles); CHK(header.vertices_count == nvertices); CHK(header.is_infinite == (i == 0)); + if(i == 0) ext_v = header.volume; + else CHK(header.volume == -ext_v); FOR_EACH(t, 0, header.primitives_count) { OK(senc3d_enclosure_get_triangle_id(enclosure, t, &gid, &side)); @@ -203,6 +206,7 @@ test(const int convention) CHK(header.unique_primitives_count == ntriangles - 1); CHK(header.vertices_count == nvertices); CHK(header.is_infinite == 1); + CHK(header.volume == 0); FOR_EACH(t, 0, header.primitives_count) { OK(senc3d_enclosure_get_triangle_id(enclosure, t, &gid, &side)); diff --git a/src/test_senc3d_inconsistant_cube.c b/src/test_senc3d_inconsistant_cube.c @@ -28,7 +28,7 @@ * 0----1' 0----1' * Front, right Back, left and * and Top faces bottom faces */ -/* Triangle #1 rotation order is not consistant with other triangles */ +/* Triangle #0 rotation order is not consistant with other triangles */ static unsigned inconsistant_box_indices[12/*#triangles*/ * 3/*#indices per triangle*/] = { 0, 1, 2, 1, 2, 3, /* Front face */ @@ -41,11 +41,21 @@ inconsistant_box_indices[12/*#triangles*/ * 3/*#indices per triangle*/] = { static const unsigned inconsistant_box_ntriangles = sizeof(inconsistant_box_indices) / (3 * sizeof(*inconsistant_box_indices)); -/* Media definitions reflect triangle #1 inconsistancy */ +/* Media definitions reflect triangle #0 inconsistancy */ static const unsigned inconsistant_medium_front[12] = { 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 }; static const unsigned inconsistant_medium_back[12] = { 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 }; +static const double inconsistant_cube_vertices[8/*#vertices*/ * 3/*#coords per vertex*/] = { + 10.0, 10.0, 10.0, + 11.0, 10.0, 10.0, + 10.0, 11.0, 10.0, + 11.0, 11.0, 10.0, + 10.0, 10.0, 11.0, + 11.0, 10.0, 11.0, + 10.0, 11.0, 11.0, + 11.0, 11.0, 11.0 +}; static void test(const int convention) @@ -68,7 +78,7 @@ test(const int convention) * but opposite sides. * What differs in this test is that triangle #0 vertices are given * in the opposite rotation order. */ - ctx.positions = box_vertices; + ctx.positions = inconsistant_cube_vertices; ctx.indices = inconsistant_box_indices; ctx.front_media = inconsistant_medium_front; ctx.back_media = inconsistant_medium_back; @@ -93,6 +103,8 @@ test(const int convention) OK(senc3d_scene_get_enclosure(scn, e, &enclosure)); OK(senc3d_enclosure_get_header(enclosure, &header)); CHK(header.enclosed_media_count == 1); + CHK(header.area == 6); + CHK(header.volume == (header.is_infinite ? -1 : 1)); OK(senc3d_enclosure_get_medium(enclosure, 0, &medium)); /* If NORMAL_FRONT the external enclosure's medium should be 0, * 1 if NORMAL_BACK */ @@ -113,7 +125,7 @@ test(const int convention) expected_side = (inconsistant_medium_front[i] == expected_medium) ? SENC3D_FRONT : SENC3D_BACK; cmp_trg(i, enclosure, - inconsistant_box_indices + (3 * i), box_vertices, + inconsistant_box_indices + (3 * i), inconsistant_cube_vertices, &same, &reversed); /* Should be made of the same triangles */ CHK(same);