stardis-solver

Solve coupled heat transfers
git clone git://git.meso-star.fr/stardis-solver.git
Log | Files | Refs | README | LICENSE

commit ae521e17da80907fb13abcc438b11400f959f6ad
parent 09c5f0f4326effabccaf1332c55acab314e9b489
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Wed, 11 Apr 2018 16:19:33 +0200

Use Star-Enc2D to analalyse 2D data

Diffstat:
Mcmake/CMakeLists.txt | 6++++--
Msrc/sdis_device_c.h | 2++
Msrc/sdis_scene.c | 795++-----------------------------------------------------------------------------
Asrc/sdis_scene_Xd.h | 767+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Msrc/sdis_scene_c.h | 13+++++++++++++
5 files changed, 799 insertions(+), 784 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -28,6 +28,7 @@ find_package(Star2D 0.1 REQUIRED) find_package(Star3D 0.4 REQUIRED) find_package(StarSP 0.7 REQUIRED) find_package(StarEnc 0.1 REQUIRED) +find_package(StarEnc2D 0.1 REQUIRED) find_package(RSys 0.6 REQUIRED) find_package(OpenMP 1.2 REQUIRED) @@ -35,7 +36,7 @@ set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${RCMAKE_SOURCE_DIR}) include(rcmake) include(rcmake_runtime) -rcmake_append_runtime_dirs(_runtime_dirs RSys Star3D StarSP StarEnc) +rcmake_append_runtime_dirs(_runtime_dirs RSys Star3D StarSP StarEnc StarEnc2D) ################################################################################ # Configure and define targets @@ -66,6 +67,7 @@ set(SDIS_FILES_INC sdis_interface_c.h sdis_medium_c.h sdis_scene_c.h + sdis_scene_Xd.h sdis_solve_Xd.h) set(SDIS_FILES_DOC COPYING README.md) @@ -80,7 +82,7 @@ add_library(sdis SHARED ${SDIS_FILES_SRC} ${SDIS_FILES_INC} ${SDIS_FILES_INC_API}) -target_link_libraries(sdis RSys Star2D Star3D StarEnc StarSP m) +target_link_libraries(sdis RSys Star2D Star3D StarEnc StarEnc2D StarSP m) if(CMAKE_COMPILER_IS_GNUCC) target_link_libraries(sdis m) endif() diff --git a/src/sdis_device_c.h b/src/sdis_device_c.h @@ -16,6 +16,8 @@ #ifndef SDIS_DEVICE_C_H #define SDIS_DEVICE_C_H +#include "sdis.h" + #include <rsys/dynamic_array.h> #include <rsys/free_list.h> #include <rsys/ref_count.h> diff --git a/src/sdis_scene.c b/src/sdis_scene.c @@ -13,791 +13,22 @@ * You should have received a copy of the GNU General Public License * along with this program. If not, see <http://www.gnu.org/licenses/>. */ -#include "sdis.h" -#include "sdis_device_c.h" -#include "sdis_interface_c.h" -#include "sdis_medium_c.h" -#include "sdis_scene_c.h" +#include "sdis_scene_Xd.h" -#include <rsys/float2.h> -#include <rsys/float3.h> -#include <rsys/double2.h> -#include <rsys/double3.h> -#include <rsys/mem_allocator.h> +/* Generate the 2D functions of the scene */ +#define SDIS_SCENE_DIMENSION 2 +#include "sdis_scene_Xd.h" -#include <senc.h> -#include <star/s2d.h> -#include <star/s3d.h> +/* Generate the 3D functions of the scene */ +#define SDIS_SCENE_DIMENSION 3 +#include "sdis_scene_Xd.h" -#include <limits.h> - -/* Context used to wrap the user geometry and interfaces to Star-Enc */ -struct geometry { - void (*indices)(const size_t iprim, size_t ids[], void*); - void (*interf)(const size_t iprim, struct sdis_interface**, void*); - void (*position)(const size_t ivert, double pos[], void*); - void* data; -}; +#include "sdis.h" +#include "sdis_scene_c.h" /******************************************************************************* * Helper function ******************************************************************************/ -/* Check that `hit' roughly lies on a vertex. For segments, a simple but - * approximative way is to test that its position have at least one barycentric - * coordinate roughly equal to 0 or 1. */ -static FINLINE int -hit_on_vertex(const struct s2d_hit* hit) -{ - const float on_vertex_eps = 1.e-4f; - float v; - ASSERT(hit && !S2D_HIT_NONE(hit)); - v = 1.f - hit->u; - return eq_epsf(hit->u, 0.f, on_vertex_eps) - || eq_epsf(hit->u, 1.f, on_vertex_eps) - || eq_epsf(v, 0.f, on_vertex_eps) - || eq_epsf(v, 1.f, on_vertex_eps); -} - -/* Check that `hit' roughly lies on an edge. For triangular primitives, a - * simple but approximative way is to test that its position have at least one - * barycentric coordinate roughly equal to 0 or 1. */ -static FINLINE int -hit_on_edge(const struct s3d_hit* hit) -{ - const float on_edge_eps = 1.e-4f; - float w; - ASSERT(hit && !S3D_HIT_NONE(hit)); - w = 1.f - hit->uv[0] - hit->uv[1]; - return eq_epsf(hit->uv[0], 0.f, on_edge_eps) - || eq_epsf(hit->uv[0], 1.f, on_edge_eps) - || eq_epsf(hit->uv[1], 0.f, on_edge_eps) - || eq_epsf(hit->uv[1], 1.f, on_edge_eps) - || eq_epsf(w, 0.f, on_edge_eps) - || eq_epsf(w, 1.f, on_edge_eps); -} - -static int -hit_filter_function_3d - (const struct s3d_hit* hit, - const float org[3], - const float dir[3], - void* ray_data, - void* filter_data) -{ - const struct s3d_hit* hit_from = ray_data; - (void)org, (void)dir, (void)filter_data; - - if(!hit_from || S3D_HIT_NONE(hit_from)) return 0; /* No filtering */ - - if(S3D_PRIMITIVE_EQ(&hit_from->prim, &hit->prim)) return 1; - - if(eq_epsf(hit->distance, 0, 1.e-6f)) { - /* If the targeted point is near of the origin, check that it lies on an - * edge shared by the 2 primitives. */ - return hit_on_edge(hit_from) && hit_on_edge(hit); - } - - return 0; -} - -static int -hit_filter_function_2d - (const struct s2d_hit* hit, - const float org[2], - const float dir[2], - void* ray_data, - void* filter_data) -{ - const struct s2d_hit* hit_from = ray_data; - (void)org, (void)dir, (void)filter_data; - - if(!hit_from || S2D_HIT_NONE(hit_from)) return 0; /* No filtering */ - - if(S2D_PRIMITIVE_EQ(&hit_from->prim, &hit->prim)) return 1; - - if(eq_epsf(hit->distance, 0, 1.e-6f)) { - /* If the targeted point is near of the origin, check that it lies on a - * vertex shared by the 2 segments. */ - return hit_on_vertex(hit_from) && hit_on_vertex(hit); - } - - return 0; -} - -static void -geometry_indices_3d(const unsigned itri, unsigned out_ids[3], void* data) -{ - struct geometry* ctx = data; - size_t ids[3]; - ASSERT(ctx && out_ids); - ctx->indices(itri, ids, ctx->data); - out_ids[0] = (unsigned)ids[0]; - out_ids[1] = (unsigned)ids[1]; - out_ids[2] = (unsigned)ids[2]; -} - -static void -geometry_media(const unsigned itri, unsigned media[2], void* data) -{ - struct geometry* ctx = data; - struct sdis_interface* interf; - ASSERT(ctx && media); - ctx->interf(itri, &interf, ctx->data); - /* FIXME check that the order in which media are returned is right */ - media[0] = medium_get_id(interf->medium_back); - media[1] = medium_get_id(interf->medium_front); -} - -static void -geometry_position_3d(const unsigned ivert, double out_pos[3], void* data) -{ - struct geometry* ctx = data; - double pos[3]; - ASSERT(ctx && out_pos); - ctx->position(ivert, pos, ctx->data); - out_pos[0] = pos[0]; - out_pos[1] = pos[1]; - out_pos[2] = pos[2]; -} - -static void -descriptor_indices_3d(const unsigned itri, unsigned ids[3], void* data) -{ - struct senc_descriptor* desc = data; - SENC(descriptor_get_global_triangle(desc, itri, ids)); -} - - -static void -descriptor_position_3d(const unsigned ivert, float out_pos[3], void* data) -{ - struct senc_descriptor* desc = data; - double pos[3]; - SENC(descriptor_get_global_vertex(desc, ivert, pos)); - out_pos[0] = (float)pos[0]; - out_pos[1] = (float)pos[1]; - out_pos[2] = (float)pos[2]; -} - -static void -enclosure_indices_3d(const unsigned itri, unsigned ids[3], void* data) -{ - struct senc_enclosure* enc = data; - SENC(enclosure_get_triangle(enc, itri, ids)); -} - -static void -enclosure_position_3d(const unsigned ivert, float out_pos[3], void* data) -{ - struct senc_enclosure* enc = data; - double pos[3]; - ASSERT(out_pos); - SENC(enclosure_get_vertex(enc, ivert, pos)); - out_pos[0] = (float)pos[0]; - out_pos[1] = (float)pos[1]; - out_pos[2] = (float)pos[2]; -} - -#if 0 -static void -get_indices_2d(const unsigned iseg, unsigned out_ids[2], void* data) -{ - struct geometry_context* ctx = data; - size_t ids[2]; - ASSERT(ctx); - ctx->indices(iseg, ids, ctx->data); - out_ids[0] = (unsigned)ids[0]; - out_ids[1] = (unsigned)ids[1]; -} - -static void -get_position_2d(const unsigned ivert, float out_pos[2], void* data) -{ - struct geometry_context* ctx = data; - double pos[2]; - ASSERT(ctx); - ctx->position(ivert, pos, ctx->data); - out_pos[0] = (float)pos[0]; - out_pos[1] = (float)pos[1]; -} - -#endif - -static void -clear_properties(struct sdis_scene* scn) -{ - size_t i; - ASSERT(scn); - FOR_EACH(i, 0, darray_interf_size_get(&scn->interfaces)) { - if(darray_interf_cdata_get(&scn->interfaces)[i]) { - SDIS(interface_ref_put(darray_interf_data_get(&scn->interfaces)[i])); - } - } - darray_interf_clear(&scn->interfaces); - darray_medium_clear(&scn->media); - darray_prim_prop_clear(&scn->prim_props); -} - -static res_T -run_analyze_3d - (struct sdis_scene* scn, - const size_t ntris, /* #triangles */ - void (*indices)(const size_t itri, size_t ids[3], void*), - void (*interf)(const size_t itri, struct sdis_interface**, void*), - const size_t nverts, /* #vertices */ - void (*position)(const size_t ivert, double pos[3], void* ctx), - void* ctx, - struct senc_descriptor** out_desc) -{ - struct geometry geom; - struct senc_device* senc = NULL; - struct senc_scene* senc_scn = NULL; - struct senc_descriptor* desc = NULL; - size_t nmedia; - size_t itri; - res_T res = RES_OK; - ASSERT(scn && ntris && indices && interf && nverts && position && out_desc); - - res = senc_device_create(scn->dev->logger, scn->dev->allocator, - scn->dev->nthreads, scn->dev->verbose, &senc); - if(res != RES_OK) goto error; - - /* Conservatively define the number of media. - * - * FIXME The number of media is going to be remove from senc_scene_create - * profile and thus the following code should be unecessary soon. */ - nmedia = 0; - FOR_EACH(itri, 0, ntris) { - struct sdis_interface* itface; - interf(itri, &itface, ctx); - nmedia = MMAX(nmedia, medium_get_id(itface->medium_front)); - nmedia = MMAX(nmedia, medium_get_id(itface->medium_back)); - } - nmedia += 1; /* +1 to define the "number of" media and not the max id */ - - res = senc_scene_create(senc, (unsigned)nmedia, &senc_scn); - if(res != RES_OK) goto error; - - /* Setup the geometry data */ - geom.indices = indices; - geom.interf = interf; - geom.position = position; - geom.data = ctx; - res = senc_scene_add_geometry - (senc_scn, (unsigned)ntris, geometry_indices_3d, geometry_media, NULL, - (unsigned)nverts, geometry_position_3d, &geom); - if(res != RES_OK) goto error; - - /* Launch the scene analyze */ - res = senc_scene_analyze(senc_scn, &desc); - if(res != RES_OK) goto error; - -exit: - if(senc) SENC(device_ref_put(senc)); - if(senc_scn) SENC(scene_ref_put(senc_scn)); - if(out_desc) *out_desc = desc; - return res; -error: - if(desc) { - SENC(descriptor_ref_put(desc)); - desc = NULL; - } - goto exit; -} - -static res_T -register_medium(struct sdis_scene* scn, struct sdis_medium* mdm) -{ - unsigned id; - size_t nmedia; - res_T res = RES_OK; - ASSERT(scn && mdm); - - /* Check that the front medium is already registered against the scene */ - id = medium_get_id(mdm); - nmedia = darray_medium_size_get(&scn->media); - if(id >= nmedia) { - res = darray_medium_resize(&scn->media, id + 1); - if(res != RES_OK) return res; - } - if(darray_medium_cdata_get(&scn->media)[id]) { - ASSERT(darray_medium_cdata_get(&scn->media)[id] == mdm); - } else { - /* Do not take a reference onto the medium since we already take a - * reference onto at least one interface that uses it, and thus that has a - * reference onto it */ - darray_medium_data_get(&scn->media)[id] = mdm; - } - return RES_OK; -} - -static res_T -setup_properties - (struct sdis_scene* scn, - struct senc_descriptor* desc, - void (*interf)(const size_t itri, struct sdis_interface**, void*), - void* ctx) -{ - unsigned itri, ntris; - res_T res = RES_OK; - ASSERT(scn && interf); - - clear_properties(scn); - - SENC(descriptor_get_global_triangles_count(desc, &ntris)); - FOR_EACH(itri, 0, ntris) { - struct prim_prop* prim_prop; - struct sdis_interface* itface; - unsigned enclosures[2]; - unsigned itri_adjusted; /* Triangle id in user space */ - unsigned id; - size_t ninterfaces; - - /* Retrieve the triangle id in user space */ - SENC(descriptor_get_global_triangle_global_id(desc, itri, &itri_adjusted)); - - /* Fetch the enclosures that the triangle splits */ - SENC(descriptor_get_global_triangle_enclosures(desc, itri, enclosures)); - - /* Fetch the interface of the primitive */ - interf(itri, &itface, ctx); - - /* Check that the interface is already registered against the scene */ - id = interface_get_id(itface); - ninterfaces = darray_interf_size_get(&scn->interfaces); - if(id >= ninterfaces) { - res = darray_interf_resize(&scn->interfaces, id + 1); - if(res != RES_OK) goto error; - } - if(darray_interf_cdata_get(&scn->interfaces)[id]) { - ASSERT(darray_interf_cdata_get(&scn->interfaces)[id] == itface); - } else { - SDIS(interface_ref_get(itface)); - darray_interf_data_get(&scn->interfaces)[id] = itface; - } - - /* Register the interface media against the scene */ - res = register_medium(scn, itface->medium_front); - if(res != RES_OK) goto error; - res = register_medium(scn, itface->medium_back); - if(res != RES_OK) goto error; - - /* Allocate primitive properties */ - res = darray_prim_prop_resize(&scn->prim_props, itri+1); - if(res != RES_OK) goto error; - - - /* Setup primitive properties */ - prim_prop = darray_prim_prop_data_get(&scn->prim_props) + itri; - prim_prop->interf = itface; - /* FIXME Ensure that the enclosure order is right one. Actually, it seems - * that Star-Enc front facing convention is reversed wrt Stardis. Faces are - * front facing when their vertex are CCW ordered */ - prim_prop->back_enclosure = enclosures[0]; - prim_prop->front_enclosure = enclosures[1]; - } - -exit: - return res; -error: - clear_properties(scn); - goto exit; -} - -#if 0 -static res_T -setup_geometry_2d - (struct sdis_scene* scn, - const size_t nsegs, /* #segments */ - void (*indices)(const size_t itri, size_t ids[2], void*), - const size_t nverts, /* #vertices */ - void (*position)(const size_t ivert, double pos[2], void* ctx), - void* ctx) -{ - struct geometry_context context; - struct s2d_shape* s2d_msh = NULL; - struct s2d_scene* s2d_scn = NULL; - struct s2d_vertex_data vdata = S2D_VERTEX_DATA_NULL; - res_T res = RES_OK; - ASSERT(scn && nsegs && indices && nverts && position); - - /* Setup the intermediary geometry context */ - context.indices = indices; - context.position = position; - context.data = ctx; - - /* Setup the vertex data */ - vdata.usage = S2D_POSITION; - vdata.type = S2D_FLOAT2; - vdata.get = get_position_2d; - - /* Create the Star-2D geometry */ - res = s2d_scene_create(scn->dev->s2d, &s2d_scn); - if(res != RES_OK) goto error; - res = s2d_shape_create_line_segments(scn->dev->s2d, &s2d_msh); - if(res != RES_OK) goto error; - res = s2d_line_segments_set_hit_filter_function - (s2d_msh, hit_filter_function_2d, NULL); - if(res != RES_OK) goto error; - res = s2d_scene_attach_shape(s2d_scn, s2d_msh); - if(res != RES_OK) goto error; - res = s2d_line_segments_setup_indexed_vertices(s2d_msh, (unsigned)nsegs, - get_indices_2d, (unsigned)nverts, &vdata, 1, &context); - if(res != RES_OK) goto error; - res = s2d_scene_view_create(s2d_scn, S2D_SAMPLE|S2D_TRACE|S2D_GET_PRIMITIVE, - &scn->s2d_view); - if(res != RES_OK) goto error; - -exit: - if(s2d_msh) S2D(shape_ref_put(s2d_msh)); - if(s2d_scn) S2D(scene_ref_put(s2d_scn)); - return res; -error: - if(scn->s2d_view) S2D(scene_view_ref_put(scn->s2d_view)); - goto exit; -} -#endif - -static res_T -setup_scene_geometry_3d(struct sdis_scene* scn, struct senc_descriptor* desc) -{ - struct s3d_shape* s3d_msh = NULL; - struct s3d_scene* s3d_scn = NULL; - struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL; - unsigned ntris, nverts; - res_T res = RES_OK; - ASSERT(scn && desc); - - SENC(descriptor_get_global_triangles_count(desc, &ntris)); - SENC(descriptor_get_global_vertices_count(desc, &nverts)); - - /* Setup the vertex data */ - vdata.usage = S3D_POSITION; - vdata.type = S3D_FLOAT3; - vdata.get = descriptor_position_3d; - - /* Create the Star-3D geometry of the whole scene */ - #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0 - CALL(s3d_scene_create(scn->dev->s3d, &s3d_scn)); - CALL(s3d_shape_create_mesh(scn->dev->s3d, &s3d_msh)); - CALL(s3d_mesh_set_hit_filter_function(s3d_msh, hit_filter_function_3d, NULL)); - CALL(s3d_scene_attach_shape(s3d_scn, s3d_msh)); - CALL(s3d_mesh_setup_indexed_vertices(s3d_msh, ntris, descriptor_indices_3d, - nverts, &vdata, 1, desc)); - CALL(s3d_scene_view_create(s3d_scn, S3D_TRACE|S3D_GET_PRIMITIVE, - &scn->s3d_view)); - #undef CALL - -exit: - if(s3d_msh) S3D(shape_ref_put(s3d_msh)); - if(s3d_scn) S3D(scene_ref_put(s3d_scn)); - return res; -error: - if(scn->s3d_view) S3D(scene_view_ref_put(scn->s3d_view)); - goto exit; -} - -static res_T -setup_enclosure_geometry_3d(struct sdis_scene* scn, struct senc_enclosure* enc) -{ - struct s3d_scene* s3d_scn = NULL; - struct s3d_shape* s3d_msh = NULL; - struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL; - const struct enclosure_header* header; - struct enclosure enc_dummy; - struct enclosure* enc_data; - unsigned itri, ntris, nverts; - res_T res = RES_OK; - ASSERT(scn && enc); - - enclosure_init(scn->dev->allocator, &enc_dummy); - - SENC(enclosure_get_header(enc, &header)); - ntris = header->triangle_count; - nverts = header->vertices_count; - - /* Register the enclosure into the scene. Use a dummy data on their - * registration. We are going to setup the data after their registration into - * the hash table in order to avoid a costly copy. In other words, the - * following hash table registration can be seen as an allocation of the - * enclosure data that are then setup. */ - res = htable_enclosure_set(&scn->enclosures, &header->enclosure_id, &enc_dummy); - if(res != RES_OK) goto error; - - /* Fetch the data of the registered enclosure */ - enc_data = htable_enclosure_find(&scn->enclosures, &header->enclosure_id); - ASSERT(enc_data != NULL); - - /* Setup the vertex data */ - vdata.usage = S3D_POSITION; - vdata.type = S3D_FLOAT3; - vdata.get = enclosure_position_3d; - - /* Create the Star-3D geometry */ - #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0 - CALL(s3d_scene_create(scn->dev->s3d, &s3d_scn)); - CALL(s3d_shape_create_mesh(scn->dev->s3d, &s3d_msh)); - CALL(s3d_scene_attach_shape(s3d_scn, s3d_msh)); - CALL(s3d_mesh_setup_indexed_vertices(s3d_msh, ntris, enclosure_indices_3d, - nverts, &vdata, 1, enc)); - CALL(s3d_scene_view_create(s3d_scn, S3D_SAMPLE, &enc_data->s3d_view)); - #undef CALL - - /* Define the identifier of the enclosure primitives in the whole scene */ - res = darray_uint_resize(&enc_data->local2global, ntris); - if(res != RES_OK) goto error; - FOR_EACH(itri, 0, ntris) { - SENC(enclosure_get_triangle_global_id - (enc, itri, darray_uint_data_get(&enc_data->local2global)+itri)); - } - -exit: - enclosure_release(&enc_dummy); - if(s3d_msh) S3D(shape_ref_put(s3d_msh)); - if(s3d_scn) S3D(scene_ref_put(s3d_scn)); - return res; -error: - htable_enclosure_erase(&scn->enclosures, &header->enclosure_id); - goto exit; -} - -static res_T -setup_enclosures_3d(struct sdis_scene* scn, struct senc_descriptor* desc) -{ - struct senc_enclosure* enc = NULL; - unsigned ienc, nencs; - res_T res = RES_OK; - ASSERT(scn && desc); - - SENC(descriptor_get_enclosure_count(desc, &nencs)); - FOR_EACH(ienc, 0, nencs) { - const struct enclosure_header* header; - const struct sdis_medium* mdm; - - SENC(descriptor_get_enclosure(desc, ienc, &enc)); - SENC(enclosure_get_header(enc, &header)); - - ASSERT(header->enclosed_medium < darray_medium_size_get(&scn->media)); - mdm = darray_medium_cdata_get(&scn->media)[header->enclosed_medium]; - ASSERT(mdm); - - /* Silently discard the solid enclosures */ - if(mdm->type == SDIS_MEDIUM_FLUID) { - res = setup_enclosure_geometry_3d(scn, enc); - if(res != RES_OK) goto error; - } - SENC(enclosure_ref_put(enc)); - enc = NULL; - } - -exit: - if(enc) SENC(enclosure_ref_put(enc)); - return res; -error: - goto exit; -} - -static res_T -scene_create - (struct sdis_device* dev, - const int is_2d, - const size_t nprims, /* #primitives */ - void (*indices)(const size_t iprim, size_t ids[], void*), - void (*interf)(const size_t iprim, struct sdis_interface** bound, void*), - const size_t nverts, /* #vertices */ - void (*position)(const size_t ivert, double pos[], void* ctx), - void* ctx, - struct sdis_scene** out_scn) -{ - struct senc_descriptor* desc = NULL; - struct sdis_scene* scn = NULL; - res_T res = RES_OK; - - if(!dev || !out_scn || !nprims || !indices || !interf || !nverts - || !position || nprims > UINT_MAX || nverts > UINT_MAX) { - res = RES_BAD_ARG; - goto error; - } - - scn = MEM_CALLOC(dev->allocator, 1, sizeof(struct sdis_scene)); - if(!scn) { - log_err(dev, "%s: could not allocate the Stardis scene.\n", FUNC_NAME); - res = RES_MEM_ERR; - goto error; - } - ref_init(&scn->ref); - SDIS(device_ref_get(dev)); - scn->dev = dev; - scn->ambient_radiative_temperature = -1; - darray_interf_init(dev->allocator, &scn->interfaces); - darray_medium_init(dev->allocator, &scn->media); - darray_prim_prop_init(dev->allocator, &scn->prim_props); - htable_enclosure_init(dev->allocator, &scn->enclosures); - - if(is_2d) FATAL("2D is not supported yet.\n"); - - res = run_analyze_3d(scn, nprims, indices, interf, nverts, position, ctx, &desc); - if(res != RES_OK) { - log_err(dev, "%s: error during the scene analysis.\n", FUNC_NAME); - goto error; - } - res = setup_properties(scn, desc, interf, ctx); - if(res != RES_OK) { - log_err(dev, "%s: could not setup the scene interfaces and their media.\n", - FUNC_NAME); - goto error; - } - res = setup_scene_geometry_3d(scn, desc); - if(res != RES_OK) { - log_err(dev, "%s: could not setup the scene geometry.\n", FUNC_NAME); - goto error; - } - res = setup_enclosures_3d(scn, desc); - if(res != RES_OK) { - log_err(dev, "%s: could not setup the enclosures.\n", FUNC_NAME); - goto error; - } - -exit: - if(out_scn) *out_scn = scn; - if(desc) SENC(descriptor_ref_put(desc)); - return res; -error: - if(scn) { - SDIS(scene_ref_put(scn)); - scn = NULL; - } - goto exit; -} - -static INLINE res_T -scene_get_medium_2d - (const struct sdis_scene* scn, - const double pos[2], - const struct sdis_medium** out_medium) -{ - const struct sdis_medium* medium = NULL; - size_t iprim, nprims; - size_t nfailures = 0; - const size_t max_failures = 10; - res_T res = RES_OK; - ASSERT(scn && pos); - - S2D(scene_view_primitives_count(scn->s2d_view, &nprims)); - FOR_EACH(iprim, 0, nprims) { - struct s2d_hit hit; - struct s2d_attrib attr; - struct s2d_primitive prim; - float s; - const float range[2] = {0.f, FLT_MAX}; - float N[2], P[2], dir[2], cos_N_dir; - s = 1.f / 3.f; - - /* Retrieve a position onto the primitive */ - S2D(scene_view_get_primitive(scn->s2d_view, (unsigned)iprim, &prim)); - S2D(primitive_get_attrib(&prim, S2D_POSITION, s, &attr)); - - /* Trace a ray from the randomw walk vertex toward the retrieved primitive - * position */ - f2_normalize(dir, f2_sub(dir, attr.value, f2_set_d2(P, pos))); - S2D(scene_view_trace_ray(scn->s2d_view, P, dir, range, NULL, &hit)); - - f2_normalize(N, hit.normal); - cos_N_dir = f2_dot(N, dir); - - /* Unforeseen error. One has to intersect a primitive ! */ - if(S2D_HIT_NONE(&hit)) { - ++nfailures; - if(nfailures < max_failures) { - continue; - } else { - res = RES_BAD_ARG; - goto error; - } - } - - if(absf(cos_N_dir) > 1.e-1f) { /* Not roughly orthognonal */ - const struct sdis_interface* interf; - interf = scene_get_interface(scn, hit.prim.prim_id); - medium = interface_get_medium - (interf, cos_N_dir < 0 ? SDIS_FRONT : SDIS_BACK); - break; - } - } - -exit: - *out_medium = medium; - return res; -error: - log_err(scn->dev, "%s: could not retrieve the medium at {%g, %g}.\n", - FUNC_NAME, SPLIT2(pos)); - goto exit; -} - -static INLINE res_T -scene_get_medium_3d - (const struct sdis_scene* scn, - const double pos[3], - const struct sdis_medium** out_medium) -{ - const struct sdis_medium* medium = NULL; - size_t iprim, nprims; - size_t nfailures = 0; - const size_t max_failures = 10; - res_T res = RES_OK; - ASSERT(scn && pos); - - S3D(scene_view_primitives_count(scn->s3d_view, &nprims)); - FOR_EACH(iprim, 0, nprims) { - struct s3d_hit hit; - struct s3d_attrib attr; - struct s3d_primitive prim; - float st[2]; - const float range[2] = {0.f, FLT_MAX}; - float N[3], P[3], dir[3], cos_N_dir; - st[0] = st[1] = 1.f / 3.f; /* Or MSVC will issue a warning */ - - /* Retrieve a position onto the primitive */ - S3D(scene_view_get_primitive(scn->s3d_view, (unsigned)iprim, &prim)); - S3D(primitive_get_attrib(&prim, S3D_POSITION, st, &attr)); - - /* Trace a ray from the randomw walk vertex toward the retrieved primitive - * position */ - f3_normalize(dir, f3_sub(dir, attr.value, f3_set_d3(P, pos))); - S3D(scene_view_trace_ray(scn->s3d_view, P, dir, range, NULL, &hit)); - - f3_normalize(N, hit.normal); - cos_N_dir = f3_dot(N, dir); - - /* Unforeseen error. One has to intersect a primitive ! */ - if(S3D_HIT_NONE(&hit)) { - ++nfailures; - if(nfailures < max_failures) { - continue; - } else { - res = RES_BAD_ARG; - goto error; - } - } - - if(absf(cos_N_dir) > 1.e-1f) { /* Not roughly orthognonal */ - const struct sdis_interface* interf; - interf = scene_get_interface(scn, hit.prim.prim_id); - medium = interface_get_medium - (interf, cos_N_dir < 0 ? SDIS_FRONT : SDIS_BACK); - break; - } - } - -exit: - *out_medium = medium; - return res; -error: - log_err(scn->dev, "%s: could not retrieve the medium at {%g, %g, %g}.\n", - FUNC_NAME, SPLIT3(pos)); - goto exit; -} - static void scene_release(ref_T * ref) { @@ -831,8 +62,8 @@ sdis_scene_create void* ctx, struct sdis_scene** out_scn) { - return scene_create - (dev, 0/*is_2D*/, ntris, indices, interf, nverts, position, ctx, out_scn); + return scene_create_3d + (dev, ntris, indices, interf, nverts, position, ctx, out_scn); } res_T @@ -846,8 +77,8 @@ sdis_scene_2d_create void* ctx, struct sdis_scene** out_scn) { - return scene_create - (dev, 1/*is_2D*/, nsegs, indices, interf, nverts, position, ctx, out_scn); + return scene_create_2d + (dev, nsegs, indices, interf, nverts, position, ctx, out_scn); } res_T diff --git a/src/sdis_scene_Xd.h b/src/sdis_scene_Xd.h @@ -0,0 +1,767 @@ +/* Copyright (C) 2016-2018 |Meso|Star> (contact@meso-star.com) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. */ + +#ifndef SDIS_SCENE_DIMENSION +#ifndef SDIS_SCENE_XD_H +#define SDIS_SCENE_XD_H + +#include "sdis_interface_c.h" +#include "sdis_medium_c.h" +#include "sdis_scene_c.h" + +#include <rsys/rsys.h> + +/* Context used to wrap the user geometry and interfaces to Star-Enc */ +struct geometry { + void (*indices)(const size_t iprim, size_t ids[], void*); + void (*interf)(const size_t iprim, struct sdis_interface**, void*); + void (*position)(const size_t ivert, double pos[], void*); + void* data; +}; + +static void +geometry_media(const unsigned iprim, unsigned media[2], void* data) +{ + struct geometry* ctx = data; + struct sdis_interface* interf; + ASSERT(ctx && media); + ctx->interf(iprim, &interf, ctx->data); + /* FIXME check that the order in which media are returned is good */ + media[0] = medium_get_id(interf->medium_back); + media[1] = medium_get_id(interf->medium_front); +} + +static res_T +register_medium(struct sdis_scene* scn, struct sdis_medium* mdm) +{ + unsigned id; + size_t nmedia; + res_T res = RES_OK; + ASSERT(scn && mdm); + + /* Check that the front medium is already registered against the scene */ + id = medium_get_id(mdm); + nmedia = darray_medium_size_get(&scn->media); + if(id >= nmedia) { + res = darray_medium_resize(&scn->media, id + 1); + if(res != RES_OK) return res; + } + if(darray_medium_cdata_get(&scn->media)[id]) { + ASSERT(darray_medium_cdata_get(&scn->media)[id] == mdm); + } else { + /* Do not take a reference onto the medium since we already take a + * reference onto at least one interface that uses it, and thus that has a + * reference onto it */ + darray_medium_data_get(&scn->media)[id] = mdm; + } + return RES_OK; +} + +static void +clear_properties(struct sdis_scene* scn) +{ + size_t i; + ASSERT(scn); + FOR_EACH(i, 0, darray_interf_size_get(&scn->interfaces)) { + if(darray_interf_cdata_get(&scn->interfaces)[i]) { + SDIS(interface_ref_put(darray_interf_data_get(&scn->interfaces)[i])); + } + } + darray_interf_clear(&scn->interfaces); + darray_medium_clear(&scn->media); + darray_prim_prop_clear(&scn->prim_props); +} + +#endif /* SDIS_SCENE_XD_H */ +#else + +#include "sdis_device_c.h" + +#include <rsys/float2.h> +#include <rsys/float3.h> +#include <rsys/double2.h> +#include <rsys/double3.h> +#include <rsys/mem_allocator.h> + +#include <limits.h> + +#if (SDIS_SCENE_DIMENSION == 2) + #include <senc2d.h> + #include <star/s2d.h> + + #define sencXd(Name) CONCAT(senc2d_, Name) + #define SENCXD SENC2D + +#elif (SDIS_SCENE_DIMENSION == 3) + #include <senc.h> + #include <star/s3d.h> + + #define sencXd(Name) CONCAT(senc_, Name) + #define SENCXD SENC +#else + #error "Invalid SDIS_SCENE_DIMENSION value." +#endif + +/* Syntactic sugar */ +#define DIM SDIS_SCENE_DIMENSION + +/* Star-XD macros generic to SDIS_SCENE_DIMENSION */ +#define sXd(Name) CONCAT(CONCAT(CONCAT(s, DIM), d_), Name) +#define SXD CONCAT(CONCAT(S, DIM), D) +#define SXD_VERTEX_DATA_NULL CONCAT(CONCAT(S,DIM),D_VERTEX_DATA_NULL) +#define SXD_POSITION CONCAT(CONCAT(S, DIM), D_POSITION) +#define SXD_TRACE CONCAT(CONCAT(S,DIM), D_TRACE) +#define SXD_SAMPLE CONCAT(CONCAT(S,DIM), D_SAMPLE) +#define SXD_GET_PRIMITIVE CONCAT(CONCAT(S,DIM), D_GET_PRIMITIVE) +#define SXD_HIT_NONE CONCAT(CONCAT(S,DIM), D_HIT_NONE) +#define SXD_PRIMITIVE_EQ CONCAT(CONCAT(S,DIM), D_PRIMITIVE_EQ) + +/* Vector macros generic to SDIS_SCENE_DIMENSION */ +#define fX(Func) CONCAT(CONCAT(CONCAT(f, DIM), _), Func) +#define fX_set_dX CONCAT(CONCAT(CONCAT(f, DIM), _set_d), DIM) + +/* Macro making generic its subimitted name to SDIS_SCENE_DIMENSION */ +#define XD(Name) CONCAT(CONCAT(CONCAT(Name, _), DIM), d) + +#if DIM == 2 +/* Check that `hit' roughly lies on a vertex. For segments, a simple but + * approximative way is to test that its position have at least one barycentric + * coordinate roughly equal to 0 or 1. */ +static FINLINE int +hit_on_vertex(const struct s2d_hit* hit) +{ + const float on_vertex_eps = 1.e-4f; + float v; + ASSERT(hit && !S2D_HIT_NONE(hit)); + v = 1.f - hit->u; + return eq_epsf(hit->u, 0.f, on_vertex_eps) + || eq_epsf(hit->u, 1.f, on_vertex_eps) + || eq_epsf(v, 0.f, on_vertex_eps) + || eq_epsf(v, 1.f, on_vertex_eps); +} + +#else /* DIM == 3 */ + +/* Check that `hit' roughly lies on an edge. For triangular primitives, a + * simple but approximative way is to test that its position have at least one + * barycentric coordinate roughly equal to 0 or 1. */ +static FINLINE int +hit_on_edge(const struct s3d_hit* hit) +{ + const float on_edge_eps = 1.e-4f; + float w; + ASSERT(hit && !S3D_HIT_NONE(hit)); + w = 1.f - hit->uv[0] - hit->uv[1]; + return eq_epsf(hit->uv[0], 0.f, on_edge_eps) + || eq_epsf(hit->uv[0], 1.f, on_edge_eps) + || eq_epsf(hit->uv[1], 0.f, on_edge_eps) + || eq_epsf(hit->uv[1], 1.f, on_edge_eps) + || eq_epsf(w, 0.f, on_edge_eps) + || eq_epsf(w, 1.f, on_edge_eps); +} +#endif + +static int +XD(hit_filter_function) + (const struct sXd(hit)* hit, + const float org[DIM], + const float dir[DIM], + void* ray_data, + void* filter_data) +{ + const struct sXd(hit)* hit_from = ray_data; + (void)org, (void)dir, (void)filter_data; + + if(!hit_from || SXD_HIT_NONE(hit_from)) return 0; /* No filtering */ + + if(SXD_PRIMITIVE_EQ(&hit_from->prim, &hit->prim)) return 1; + + if(eq_epsf(hit->distance, 0, 1.e-6f)) { + /* If the targeted point is near of the origin, check that it lies on an + * edge/vertex shared by the 2 primitives. */ +#if DIM == 2 + return hit_on_vertex(hit_from) && hit_on_vertex(hit); +#else + return hit_on_edge(hit_from) && hit_on_edge(hit); +#endif + } + + return 0; +} + +static void +XD(geometry_indices)(const unsigned iprim, unsigned out_ids[DIM], void* data) +{ + struct geometry* ctx = data; + size_t ids[DIM]; + int i; + ASSERT(ctx && out_ids); + ctx->indices(iprim, ids, ctx->data); + FOR_EACH(i, 0, DIM) out_ids[i] = (unsigned)ids[i]; +} + +static void +XD(geometry_position)(const unsigned ivert, double out_pos[DIM], void* data) +{ + struct geometry* ctx = data; + double pos[DIM]; + int i; + ASSERT(ctx && out_pos); + ctx->position(ivert, pos, ctx->data); + FOR_EACH(i, 0, DIM) out_pos[i] = pos[i]; +} + +static void +XD(descriptor_indices)(const unsigned iprim, unsigned ids[DIM], void* data) +{ + struct sencXd(descriptor)* desc = data; +#if DIM == 2 + SENCXD(descriptor_get_global_segment(desc, iprim, ids)); +#else + SENCXD(descriptor_get_global_triangle(desc, iprim, ids)); +#endif +} + +static void +XD(descriptor_position)(const unsigned ivert, float out_pos[DIM], void* data) +{ + struct sencXd(descriptor)* desc = data; + double pos[3]; + int i; + SENCXD(descriptor_get_global_vertex(desc, ivert, pos)); + FOR_EACH(i, 0, DIM) out_pos[i] = (float)pos[i]; +} + +static void +XD(enclosure_indices)(const unsigned iprim, unsigned ids[DIM], void* data) +{ + struct sencXd(enclosure)* enc = data; +#if DIM == 2 + SENCXD(enclosure_get_segment(enc, iprim, ids)); +#else + SENCXD(enclosure_get_triangle(enc, iprim, ids)); +#endif +} + +static void +XD(enclosure_position)(const unsigned ivert, float out_pos[DIM], void* data) +{ + struct sencXd(enclosure)* enc = data; + double pos[DIM]; + int i; + ASSERT(out_pos); + SENCXD(enclosure_get_vertex(enc, ivert, pos)); + FOR_EACH(i, 0, DIM) out_pos[i] = (float)pos[i]; +} + +static res_T +XD(run_analyze) + (struct sdis_scene* scn, + const size_t nprims, /* #primitives */ + void (*indices)(const size_t iprim, size_t ids[], void*), + void (interf)(const size_t iprim, struct sdis_interface**, void*), + const size_t nverts, /* #vertices */ + void (*position)(const size_t ivert, double pos[], void*), + void* ctx, + struct sencXd(descriptor)** out_desc) +{ + struct geometry geom; + struct sencXd(device)* senc = NULL; + struct sencXd(scene)* senc_scn = NULL; + struct sencXd(descriptor)* desc = NULL; + size_t nmedia; + size_t iprim; + res_T res = RES_OK; + ASSERT(scn && nprims && indices && interf && nverts && position && out_desc); + + res = sencXd(device_create)(scn->dev->logger, scn->dev->allocator, + scn->dev->nthreads, scn->dev->verbose, &senc); + if(res != RES_OK) goto error; + + /* Conservatively define the number of media. + * + * FIXME The number of media is going to be remove from senc_scene_create + * profile and thus the following code should be unecessary soon. */ + nmedia = 0; + FOR_EACH(iprim, 0, nprims) { + struct sdis_interface* itface; + interf(iprim, &itface, ctx); + nmedia = MMAX(nmedia, medium_get_id(itface->medium_front)); + nmedia = MMAX(nmedia, medium_get_id(itface->medium_back)); + } + nmedia += 1; /* +1 to define the "number of" media and not the max id */ + + res = sencXd(scene_create)(senc, (unsigned)nmedia, &senc_scn); + if(res != RES_OK) goto error; + + /* Setup the geometry data */ + geom.indices = indices; + geom.interf = interf; + geom.position = position; + geom.data = ctx; + res = sencXd(scene_add_geometry) + (senc_scn, (unsigned)nprims, XD(geometry_indices), geometry_media, NULL, + (unsigned)nverts, XD(geometry_position), &geom); + if(res != RES_OK) goto error; + + /* Launch the scene analyze */ + res = sencXd(scene_analyze)(senc_scn, &desc); + if(res != RES_OK) goto error; + +exit: + if(senc) SENCXD(device_ref_put(senc)); + if(senc_scn) SENCXD(scene_ref_put(senc_scn)); + if(out_desc) *out_desc = desc; + return res; +error: + if(desc) { + SENCXD(descriptor_ref_put(desc)); + desc = NULL; + } + goto exit; + +} + +static res_T +XD(setup_properties) + (struct sdis_scene* scn, + struct sencXd(descriptor)* desc, + void (*interf)(const size_t itri, struct sdis_interface**, void*), + void* ctx) +{ + unsigned iprim, nprims; + res_T res = RES_OK; + ASSERT(scn && interf); + + clear_properties(scn); + +#if DIM == 2 + SENCXD(descriptor_get_global_segments_count(desc, &nprims)); +#else + SENCXD(descriptor_get_global_triangles_count(desc, &nprims)); +#endif + FOR_EACH(iprim, 0, nprims) { + struct prim_prop* prim_prop; + struct sdis_interface* itface; + unsigned enclosures[2]; + unsigned iprim_adjusted; /* Primitive id in user space */ + unsigned id; + size_t ninterfaces; + +#if DIM == 2 + /* Retrieve the triangle id in user space */ + SENCXD(descriptor_get_global_segment_global_id(desc, iprim, &iprim_adjusted)); + /* Fetch the enclosures that the segment splits */ + SENCXD(descriptor_get_global_segment_enclosures(desc, iprim, enclosures)); +#else + /* Retrieve the triangle id in user space */ + SENCXD(descriptor_get_global_triangle_global_id(desc, iprim, &iprim_adjusted)); + /* Fetch the enclosures that the triangle splits */ + SENCXD(descriptor_get_global_triangle_enclosures(desc, iprim, enclosures)); +#endif + + /* Fetch the interface of the primitive */ + interf(iprim, &itface, ctx); + + /* Check that the interface is already registered against the scene */ + id = interface_get_id(itface); + ninterfaces = darray_interf_size_get(&scn->interfaces); + if(id >= ninterfaces) { + res = darray_interf_resize(&scn->interfaces, id + 1); + if(res != RES_OK) goto error; + } + if(darray_interf_cdata_get(&scn->interfaces)[id]) { + ASSERT(darray_interf_cdata_get(&scn->interfaces)[id] == itface); + } else { + SDIS(interface_ref_get(itface)); + darray_interf_data_get(&scn->interfaces)[id] = itface; + } + + /* Register the interface media against the scene */ + res = register_medium(scn, itface->medium_front); + if(res != RES_OK) goto error; + res = register_medium(scn, itface->medium_back); + if(res != RES_OK) goto error; + + /* Allocate primitive properties */ + res = darray_prim_prop_resize(&scn->prim_props, iprim+1); + if(res != RES_OK) goto error; + + /* Setup primitive properties */ + prim_prop = darray_prim_prop_data_get(&scn->prim_props) + iprim; + prim_prop->interf = itface; + /* FIXME Ensure that the enclosure order is right one. Actually, it seems + * that Star-Enc front facing convention is reversed wrt Stardis. Faces are + * front facing when their vertex are CCW ordered */ + prim_prop->back_enclosure = enclosures[0]; + prim_prop->front_enclosure = enclosures[1]; + } + +exit: + return res; +error: + clear_properties(scn); + goto exit; +} + +static res_T +XD(setup_scene_geometry)(struct sdis_scene* scn, struct sencXd(descriptor)* desc) +{ + struct sXd(device)* sXd_dev = NULL; + struct sXd(shape)* sXd_shape = NULL; + struct sXd(scene)* sXd_scn = NULL; + struct sXd(vertex_data) vdata = SXD_VERTEX_DATA_NULL; + unsigned nprims, nverts; + res_T res = RES_OK; + ASSERT(scn && desc); + + SENCXD(descriptor_get_global_vertices_count(desc, &nverts)); + + /* Setup the vertex data */ + vdata.usage = SXD_POSITION; +#if DIM == 2 + vdata.type = S2D_FLOAT2; +#else + vdata.type = S3D_FLOAT3; +#endif + vdata.get = XD(descriptor_position); + + /* Create the Star-XD geometry of the whole scene */ + #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0 +#if DIM == 2 + sXd_dev = scn->dev->s2d; + SENCXD(descriptor_get_global_segments_count(desc, &nprims)); + CALL(sXd(shape_create_line_segments)(sXd_dev, &sXd_shape)); + CALL(sXd(line_segments_set_hit_filter_function)(sXd_shape, + XD(hit_filter_function), NULL)); + CALL(sXd(line_segments_setup_indexed_vertices)(sXd_shape, nprims, + XD(descriptor_indices), nverts, &vdata, 1, desc)); +#else + sXd_dev = scn->dev->s3d; + SENCXD(descriptor_get_global_triangles_count(desc, &nprims)); + CALL(sXd(shape_create_mesh)(sXd_dev, &sXd_shape)); + CALL(sXd(mesh_set_hit_filter_function)(sXd_shape, XD(hit_filter_function), NULL)); + CALL(sXd(mesh_setup_indexed_vertices)(sXd_shape, nprims, XD(descriptor_indices), + nverts, &vdata, 1, desc)); +#endif + CALL(sXd(scene_create)(sXd_dev, &sXd_scn)); + CALL(sXd(scene_attach_shape)(sXd_scn, sXd_shape)); + CALL(sXd(scene_view_create)(sXd_scn, SXD_TRACE|SXD_GET_PRIMITIVE, + &scn->sXd(view))); + #undef CALL + +exit: + if(sXd_shape) SXD(shape_ref_put(sXd_shape)); + if(sXd_scn) SXD(scene_ref_put(sXd_scn)); + return res; +error: + if(scn->sXd(view)) SXD(scene_view_ref_put(scn->sXd(view))); + goto exit; +} + +static res_T +XD(setup_enclosure_geometry)(struct sdis_scene* scn, struct sencXd(enclosure)* enc) +{ + struct sXd(device)* sXd_dev = NULL; + struct sXd(scene)* sXd_scn = NULL; + struct sXd(shape)* sXd_shape = NULL; + struct sXd(vertex_data) vdata = SXD_VERTEX_DATA_NULL; + struct enclosure enc_dummy; + struct enclosure* enc_data; + unsigned iprim, nprims, nverts; +#if DIM == 2 + const struct enclosure2d_header* header; +#else + const struct enclosure_header* header; +#endif + res_T res = RES_OK; + ASSERT(scn && enc); + + enclosure_init(scn->dev->allocator, &enc_dummy); + + SENCXD(enclosure_get_header(enc, &header)); +#if DIM == 2 + sXd_dev = scn->dev->s2d; + nprims = header->segment_count; +#else + sXd_dev = scn->dev->s3d; + nprims = header->triangle_count; +#endif + nverts = header->vertices_count; + + /* Register the enclosure into the scene. Use a dummy data on their + * registration. We are going to setup the data after their registration into + * the hash table in order to avoid a costly copy. In other words, the + * following hash table registration can be seen as an allocation of the + * enclosure data that are then setup. */ + res = htable_enclosure_set(&scn->enclosures, &header->enclosure_id, &enc_dummy); + if(res != RES_OK) goto error; + + /* Fetch the data of the registered enclosure */ + enc_data = htable_enclosure_find(&scn->enclosures, &header->enclosure_id); + ASSERT(enc_data != NULL); + + /* Setup the vertex data */ + vdata.usage = SXD_POSITION; +#if DIM == 2 + vdata.type = S2D_FLOAT2; +#else + vdata.type = S2D_FLOAT3; +#endif + vdata.get = XD(enclosure_position); + + /* Create the Star-XD geometry */ + #define CALL(Func) { if(RES_OK != (res = Func)) goto error; } (void)0 +#if DIM == 2 + CALL(sXd(shape_create_line_segments)(sXd_dev, &sXd_shape)); + CALL(sXd(line_segments_setup_indexed_vertices)(sXd_shape, nprims, + XD(enclosure_indices), nverts, &vdata, 1, enc)); +#else + CALL(sXd(shape_create_mesh)(sXd_dev, &sXd_shape)); + CALL(sXd(mesh_setup_indexed_vertices)(sXd_shape, nprims, XD(enclosure_indices), + nverts, &vdata, 1, enc)); +#endif + + CALL(sXd(scene_create)(sXd_dev, &sXd_scn)); + CALL(sXd(scene_attach_shape)(sXd_scn, sXd_shape)); + CALL(sXd(scene_view_create)(sXd_scn, SXD_SAMPLE, &enc_data->sXd(view))); + #undef CALL + + /* Define the identifier of the enclosure primitives in the whole scene */ + res = darray_uint_resize(&enc_data->local2global, nprims); + if(res != RES_OK) goto error; + FOR_EACH(iprim, 0, nprims) { +#if DIM == 2 + SENCXD(enclosure_get_segment_global_id + (enc, iprim, darray_uint_data_get(&enc_data->local2global)+iprim)); +#else + SENCXD(enclosure_get_triangle_global_id + (enc, iprim, darray_uint_data_get(&enc_data->local2global)+iprim)); +#endif + } + +exit: + enclosure_release(&enc_dummy); + if(sXd_shape) SXD(shape_ref_put(sXd_shape)); + if(sXd_scn) SXD(scene_ref_put(sXd_scn)); + return res; +error: + htable_enclosure_erase(&scn->enclosures, &header->enclosure_id); + goto exit; +} + +static res_T +XD(setup_enclosures)(struct sdis_scene* scn, struct sencXd(descriptor)* desc) +{ + struct sencXd(enclosure)* enc = NULL; + unsigned ienc, nencs; + res_T res = RES_OK; + ASSERT(scn && desc); + + SENCXD(descriptor_get_enclosure_count(desc, &nencs)); + FOR_EACH(ienc, 0, nencs) { +#if DIM == 2 + const struct enclosure2d_header* header; +#else + const struct enclosure_header* header; +#endif + const struct sdis_medium* mdm; + + SENCXD(descriptor_get_enclosure(desc, ienc, &enc)); + SENCXD(enclosure_get_header(enc, &header)); + + ASSERT(header->enclosed_medium < darray_medium_size_get(&scn->media)); + mdm = darray_medium_cdata_get(&scn->media)[header->enclosed_medium]; + ASSERT(mdm); + + /* Silently discard the solid enclosures */ + if(mdm->type == SDIS_MEDIUM_FLUID) { + res = XD(setup_enclosure_geometry)(scn, enc); + if(res != RES_OK) goto error; + } + SENCXD(enclosure_ref_put(enc)); + enc = NULL; + } + +exit: + if(enc) SENCXD(enclosure_ref_put(enc)); + return res; +error: + goto exit; +} + +static res_T +XD(scene_create) + (struct sdis_device* dev, + const size_t nprims, /* #primitives */ + void (*indices)(const size_t iprim, size_t ids[], void*), + void (*interf)(const size_t iprim, struct sdis_interface** bound, void*), + const size_t nverts, /* #vertices */ + void (*position)(const size_t ivert, double pos[], void* ctx), + void* ctx, + struct sdis_scene** out_scn) +{ + struct sencXd(descriptor)* desc = NULL; + struct sdis_scene* scn = NULL; + res_T res = RES_OK; + + if(!dev || !out_scn || !nprims || !indices || !interf || !nverts + || !position || nprims > UINT_MAX || nverts > UINT_MAX) { + res = RES_BAD_ARG; + goto error; + } + + scn = MEM_CALLOC(dev->allocator, 1, sizeof(struct sdis_scene)); + if(!scn) { + log_err(dev, "%s: could not allocate the Stardis scene.\n", FUNC_NAME); + res = RES_MEM_ERR; + goto error; + } + ref_init(&scn->ref); + SDIS(device_ref_get(dev)); + scn->dev = dev; + scn->ambient_radiative_temperature = -1; + darray_interf_init(dev->allocator, &scn->interfaces); + darray_medium_init(dev->allocator, &scn->media); + darray_prim_prop_init(dev->allocator, &scn->prim_props); + htable_enclosure_init(dev->allocator, &scn->enclosures); + + res = XD(run_analyze)(scn, nprims, indices, interf, nverts, position, ctx, &desc); + if(res != RES_OK) { + log_err(dev, "%s: error during the scene analysis.\n", FUNC_NAME); + goto error; + } + res = XD(setup_properties)(scn, desc, interf, ctx); + if(res != RES_OK) { + log_err(dev, "%s: could not setup the scene interfaces and their media.\n", + FUNC_NAME); + goto error; + } + res = XD(setup_scene_geometry)(scn, desc); + if(res != RES_OK) { + log_err(dev, "%s: could not setup the scene geometry.\n", FUNC_NAME); + goto error; + } + res = XD(setup_enclosures)(scn, desc); + if(res != RES_OK) { + log_err(dev, "%s: could not setup the enclosures.\n", FUNC_NAME); + goto error; + } + +exit: + if(out_scn) *out_scn = scn; + if(desc) SENCXD(descriptor_ref_put(desc)); + return res; +error: + if(scn) { + SDIS(scene_ref_put(scn)); + scn = NULL; + } + goto exit; +} + +static INLINE res_T +XD(scene_get_medium) + (const struct sdis_scene* scn, + const double pos[3], + const struct sdis_medium** out_medium) +{ + const struct sdis_medium* medium = NULL; + size_t iprim, nprims; + size_t nfailures = 0; + const size_t max_failures = 10; + res_T res = RES_OK; + ASSERT(scn && pos); + + SXD(scene_view_primitives_count(scn->sXd(view), &nprims)); + FOR_EACH(iprim, 0, nprims) { + struct sXd(hit) hit; + struct sXd(attrib) attr; + struct sXd(primitive) prim; + const float range[2] = {0.f, FLT_MAX}; + float N[DIM], P[DIM], dir[DIM], cos_N_dir; +#if DIM == 2 + float st; + st = 1.f / 3.f; +#else + float st[2]; + st[0] = st[1] = 1.f / 3.f; /* Or MSVC will issue a warning */ +#endif + + /* Retrieve a position onto the primitive */ + SXD(scene_view_get_primitive(scn->sXd(view), (unsigned)iprim, &prim)); + SXD(primitive_get_attrib(&prim, SXD_POSITION, st, &attr)); + + /* Trace a ray from the randomw walk vertex toward the retrieved primitive + * position */ + fX(normalize)(dir, fX(sub)(dir, attr.value, fX_set_dX(P, pos))); + SXD(scene_view_trace_ray(scn->sXd(view), P, dir, range, NULL, &hit)); + + fX(normalize)(N, hit.normal); + cos_N_dir = fX(dot)(N, dir); + + /* Unforeseen error. One has to intersect a primitive ! */ + if(SXD_HIT_NONE(&hit)) { + ++nfailures; + if(nfailures < max_failures) { + continue; + } else { + res = RES_BAD_ARG; + goto error; + } + } + + if(absf(cos_N_dir) > 1.e-1f) { /* Not roughly orthognonal */ + const struct sdis_interface* interf; + interf = scene_get_interface(scn, hit.prim.prim_id); + medium = interface_get_medium + (interf, cos_N_dir < 0 ? SDIS_FRONT : SDIS_BACK); + break; + } + } + +exit: + *out_medium = medium; + return res; +error: +#if DIM == 2 + log_err(scn->dev, "%s: could not retrieve the medium at {%g, %g}.\n", + FUNC_NAME, SPLIT2(pos)); +#else + log_err(scn->dev, "%s: could not retrieve the medium at {%g, %g, %g}.\n", + FUNC_NAME, SPLIT3(pos)); +#endif + goto exit; +} + +#undef SDIS_SCENE_DIMENSION +#undef DIM +#undef sencXd +#undef SENCXD +#undef sXd +#undef SXD +#undef SXD_VERTEX_DATA_NULL +#undef SXD_POSITION +#undef SXD_FLOAT3 +#undef SXD_TRACE +#undef SXD_TRACE +#undef SXD_GET_PRIMITIVE +#undef SXD_HIT_NONE +#undef SXD_PRIMITIVE_EQ +#undef fX +#undef fX_set_dX +#undef XD + +#endif /* !SDIS_SCENE_DIMENSION */ diff --git a/src/sdis_scene_c.h b/src/sdis_scene_c.h @@ -16,6 +16,7 @@ #ifndef SDIS_SCENE_C_H #define SDIS_SCENE_C_H +#include <star/s2d.h> #include <star/s3d.h> #include <rsys/dynamic_array_uint.h> @@ -54,6 +55,7 @@ medium_init(struct mem_allocator* allocator, struct sdis_medium** medium) } struct enclosure { + struct s2d_scene_view* s2d_view; struct s3d_scene_view* s3d_view; /* Map the id of the enclosure primitives to their primitive id into the * whole scene */ @@ -64,6 +66,7 @@ static INLINE void enclosure_init(struct mem_allocator* allocator, struct enclosure* enc) { ASSERT(allocator && enc); + enc->s2d_view = NULL; enc->s3d_view = NULL; darray_uint_init(allocator, &enc->local2global); } @@ -71,6 +74,7 @@ enclosure_init(struct mem_allocator* allocator, struct enclosure* enc) static INLINE void enclosure_release(struct enclosure* enc) { + if(enc->s2d_view) S2D(scene_view_ref_put(enc->s2d_view)); if(enc->s3d_view) S3D(scene_view_ref_put(enc->s3d_view)); darray_uint_release(&enc->local2global); } @@ -82,6 +86,10 @@ enclosure_copy(struct enclosure* dst, const struct enclosure* src) S3D(scene_view_ref_get(src->s3d_view)); dst->s3d_view = src->s3d_view; } + if(src->s2d_view) { + S2D(scene_view_ref_get(src->s2d_view)); + dst->s2d_view = src->s2d_view; + } return darray_uint_copy(&dst->local2global, &src->local2global); } @@ -96,6 +104,11 @@ enclosure_copy_and_release(struct enclosure* dst, struct enclosure* src) dst->s3d_view = src->s3d_view; src->s3d_view = NULL; } + if(src->s2d_view) { + /* Only transfer ownership */ + dst->s2d_view = src->s2d_view; + src->s2d_view = NULL; + } return RES_OK; }