star-3d

Surface structuring for efficient 3D geometric queries
git clone git://git.meso-star.fr/star-3d.git
Log | Files | Refs | README | LICENSE

commit 7dba474716e2377995f66f532085952246704a9b
parent cdf047d463ab7cfb6b709dc389cc59ea271fc84e
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Fri, 11 Oct 2019 15:09:20 +0200

Merge branch 'feature_point_query' into develop

Diffstat:
Mcmake/CMakeLists.txt | 5++++-
Msrc/s3d.h | 44+++++++++++++++++++++++++++++++-------------
Msrc/s3d_geometry.c | 20++++----------------
Msrc/s3d_geometry.h | 2+-
Msrc/s3d_scene_view.c | 252+------------------------------------------------------------------------------
Asrc/s3d_scene_view_closest_point.c | 440+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Asrc/s3d_scene_view_trace_ray.c | 302++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Msrc/s3d_sphere.h | 22++++++++++++++++++++++
Asrc/test_s3d_closest_point.c | 822+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
9 files changed, 1627 insertions(+), 282 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -39,7 +39,7 @@ option(NO_TEST "Disable the test" OFF) ################################################################################ # Check dependencies ################################################################################ -find_package(Embree 3 REQUIRED) +find_package(Embree 3.6 REQUIRED) find_package(RCMake 0.2.2 REQUIRED) find_package(RSys 0.6 REQUIRED) @@ -69,6 +69,8 @@ set(S3D_FILES_SRC s3d_primitive.c s3d_scene.c s3d_scene_view.c + s3d_scene_view_closest_point.c + s3d_scene_view_trace_ray.c s3d_shape.c s3d_sphere.c) set(S3D_FILES_INC_API s3d.h) @@ -131,6 +133,7 @@ if(NOT NO_TEST) endfunction() new_test(test_s3d_accel_struct_conf) + new_test(test_s3d_closest_point) new_test(test_s3d_device) new_test(test_s3d_sampler) new_test(test_s3d_sample_sphere) diff --git a/src/s3d.h b/src/s3d.h @@ -146,7 +146,7 @@ struct s3d_hit { struct s3d_primitive prim; /* Intersected primitive */ float normal[3]; /* Un-normalized geometry normal */ float uv[2]; /* Barycentric coordinates of the hit onto `prim' */ - float distance; /* Hit distance from the ray origin */ + float distance; /* Hit distance from the query origin */ }; /* Constant defining a NULL intersection. Should be used to initialize a hit */ @@ -205,15 +205,16 @@ static const struct s3d_accel_struct_conf S3D_ACCEL_STRUCT_CONF_DEFAULT = S3D_ACCEL_STRUCT_CONF_DEFAULT__; /* Filter function data type. One can define such function to discard - * intersections along a ray with respect to user defined criteria, e.g.: - * masked/transparent primitive, etc. Return 0 or the intersection is not - * discarded and a value not equal to zero otherwise. */ + * intersections along a ray or the result of a closest point query with + * respect to user defined criteria, e.g.: masked/transparent primitive, etc. + * Return 0 or the intersection is not discarded and a value not equal to zero + * otherwise. */ typedef int (*s3d_hit_filter_function_T) (const struct s3d_hit* hit, - const float ray_org[3], - const float ray_dir[3], - void* ray_data, /* User data submitted on trace ray(s) invocation */ + const float org[3], + const float dir[3], /* Direction from `org' to `hit' */ + void* query_data, /* User data submitted on query invocation */ void* filter_data); /* Data defined on the setup of the filter function */ /* Forward declaration of s3d opaque data types */ @@ -339,6 +340,11 @@ s3d_scene_view_get_mask (struct s3d_scene_view* scnview, int* mask); +/* Trace a ray into the scene and return the closest intersection along it. The + * ray is defined by `origin' + t*`direction' = 0 with t in [`range[0]', + * `range[1]'). Note that if a range is degenerated (i.e. `range[0]' >= + * `range[1]') then the ray is not traced and `hit' is set to S3D_HIT_NULL. Can + * be called only if the scnview was created with the S3D_TRACE flag. */ S3D_API res_T s3d_scene_view_trace_ray (struct s3d_scene_view* scnview, @@ -348,12 +354,8 @@ s3d_scene_view_trace_ray void* ray_data, /* User ray data sent to the hit filter func. May be NULL */ struct s3d_hit* hit); -/* Trace a bundle of rays into the scene and return the closest intersection - * along them. The rays are defined by `origin' + t*`direction' = 0 with t in - * [`range[0]', `range[1]'). Note that if a range is degenerated (i.e. - * `range[0]' >= `range[1]') then its associated ray is not traced and `hit' is - * set to S3D_HIT_NULL. Can be called only if the scnview was created with the - * S3D_TRACE flag. */ +/* Trace a bundle of rays into the scene. Can be called only if the scnview was + * created with the S3D_TRACE flag. */ S3D_API res_T s3d_scene_view_trace_rays (struct s3d_scene_view* scnview, @@ -366,6 +368,22 @@ s3d_scene_view_trace_rays const size_t sizeof_ray_data, /* Size in Bytes of *one* ray data */ struct s3d_hit* hits); +/* Return the point onto the scene surfaces that is the closest of the + * submitted `pos'. The `radius' parameter defines the maximum search distance + * around `pos'. Each candidate position are internally filtered by the + * hit_filter_function attached to the corresponding shape; the user can thus + * reject a candidate position according to its own criteria. This function can + * be called only if the scnview was created with the S3D_TRACE flag which is + * actually the flag used to tell Star-3D to internally build an acceleration + * structure on which this function relies. */ +S3D_API res_T +s3d_scene_view_closest_point + (struct s3d_scene_view* scnview, + const float pos[3], /* Position to query */ + const float radius, /* Search distance in [0, radius[ */ + void* query_data, /* User data sent to the hit filter func. May be NULL */ + struct s3d_hit* hit); + /* Uniformly sample the scene and return the sampled primitive and its sample * uv position. Can be called only if the scnview was created with the * S3D_SAMPLE flag */ diff --git a/src/s3d_geometry.c b/src/s3d_geometry.c @@ -51,9 +51,8 @@ sphere_ray_hit_setup struct RTCHitN* hitN; struct RTCHit hit; struct RTCRay ray; - float cos_theta; float Ng[3]; - float u, v; + float uv[2]; size_t i; ASSERT(args && args->primID == 0 && args->N == 1 && args->valid[0] != 0); @@ -71,24 +70,13 @@ sphere_ray_hit_setup Ng[2] = ray.dir_z*tfar + ray.org_z - geom->data.sphere->pos[2]; f3_normalize(Ng, Ng); - cos_theta = Ng[2]; - - v = (1.f - cos_theta) * 0.5f; - if(absf(cos_theta) == 1) { - u = 0; - } else if(eq_epsf(Ng[0], 0.f, 1.e-6f)) { - u = Ng[1] > 0 ? 0.25f : 0.75f; - } else { - double phi = atan2f(Ng[1], Ng[0]); /* phi in [-PI, PI] */ - if(phi < 0) phi = 2*PI + phi; /* phi in [0, 2PI] */ - u = (float)(phi / (2*PI)); - } + sphere_normal_to_uv(Ng, uv); hit.Ng_x = Ng[0]; hit.Ng_y = Ng[1]; hit.Ng_z = Ng[2]; - hit.u = u; - hit.v = v; + hit.u = uv[0]; + hit.v = uv[1]; hit.primID = 0; hit.geomID = geom->rtc_id; FOR_EACH(i, 0, RTC_MAX_INSTANCE_LEVEL_COUNT) { diff --git a/src/s3d_geometry.h b/src/s3d_geometry.h @@ -35,6 +35,7 @@ #include "s3d.h" #include "s3d_backend.h" +#include <rsys/float3.h> #include <rsys/ref_count.h> enum geometry_type { @@ -100,4 +101,3 @@ geometry_rtc_sphere_intersect (const struct RTCIntersectFunctionNArguments* args); #endif /* S3D_GEOMETRY_H */ - diff --git a/src/s3d_scene_view.c b/src/s3d_scene_view.c @@ -41,14 +41,6 @@ #include <rsys/float33.h> #include <rsys/mem_allocator.h> -struct intersect_context { - struct RTCIntersectContext rtc; - struct s3d_scene_view* scnview; - void* data; /* Per ray user defined data */ - float ws_org[3]; /* World space ray origin */ - float ws_dir[3]; /* World space ray direction */ -}; - /******************************************************************************* * Helper functions ******************************************************************************/ @@ -138,100 +130,6 @@ on_shape_detach } } -static INLINE void -hit_setup - (struct s3d_scene_view* scnview, - const struct RTCRayHit* ray_hit, - struct s3d_hit* hit) -{ - float w; - char flip_surface = 0; - - ASSERT(scnview && hit && ray_hit); - - if(ray_hit->hit.geomID == RTC_INVALID_GEOMETRY_ID) { /* No hit */ - *hit = S3D_HIT_NULL; - return; - } - - hit->normal[0] = ray_hit->hit.Ng_x; - hit->normal[1] = ray_hit->hit.Ng_y; - hit->normal[2] = ray_hit->hit.Ng_z; - hit->distance = ray_hit->ray.tfar; - - if(ray_hit->hit.instID[0] == RTC_INVALID_GEOMETRY_ID) { - struct geometry* geom_shape; - geom_shape = scene_view_geometry_from_embree_id(scnview, ray_hit->hit.geomID); - hit->prim.shape__ = geom_shape; - hit->prim.inst__ = NULL; - hit->prim.prim_id = ray_hit->hit.primID; - hit->prim.geom_id = geom_shape->name; - hit->prim.inst_id = S3D_INVALID_ID; - hit->prim.scene_prim_id = /* Compute the "scene space" primitive id */ - hit->prim.prim_id /* Mesh space */ - + geom_shape->scene_prim_id_offset; /* Scene space */ - - } else { /* The hit shape is instantiated */ - /* Retrieve the hit instance */ - struct geometry* geom_inst; - struct geometry* geom_shape; - float transform[9]; - geom_inst = scene_view_geometry_from_embree_id - (scnview, ray_hit->hit.instID[0]); - geom_shape = scene_view_geometry_from_embree_id - (geom_inst->data.instance->scnview, ray_hit->hit.geomID); - hit->prim.shape__ = geom_shape; - hit->prim.inst__ = geom_inst; - hit->prim.prim_id = ray_hit->hit.primID; - hit->prim.geom_id = geom_shape->name; - hit->prim.inst_id = geom_inst->name; - hit->prim.scene_prim_id = /* Compute the "scene space" primitive id */ - hit->prim.prim_id /* Shape space */ - + geom_shape->scene_prim_id_offset /* Inst space */ - + geom_inst->scene_prim_id_offset; /* Scene space */ - - flip_surface = geom_inst->flip_surface; - ASSERT(hit->prim.inst__); - ASSERT(((struct geometry*)hit->prim.inst__)->type == GEOM_INSTANCE); - - /* Transform the normal in world space */ - f33_invtrans(transform, geom_inst->data.instance->transform); - f33_mulf3(hit->normal, transform, hit->normal); - } - ASSERT(hit->prim.shape__); - ASSERT(((struct geometry*)hit->prim.shape__)->type == GEOM_MESH - ||((struct geometry*)hit->prim.shape__)->type == GEOM_SPHERE); - - hit->uv[0] = ray_hit->hit.u; - hit->uv[1] = ray_hit->hit.v; - - if(((struct geometry*)hit->prim.shape__)->type == GEOM_MESH) { - w = 1.f - hit->uv[0] - hit->uv[1]; - ASSERT(w <= 1.f); /* This may not occurs */ - if(w < 0.f) { /* Handle precision error */ - if(hit->uv[0] > hit->uv[1]) hit->uv[0] += w; - else hit->uv[1] += w; - w = 0.f; - } - - /* Embree stores on the u and v ray parameters the barycentric coordinates of - * the hit with respect to the second and third triangle vertices, - * respectively. The following code computes the barycentric coordinates of - * the hit for the first and second triangle vertices */ - hit->uv[1] = hit->uv[0]; - hit->uv[0] = w; - - /* In Embree3 the normal orientation is flipped wrt to Star-3D convention */ - #if RTC_VERSION_MAJOR >= 3 - f3_minus(hit->normal, hit->normal); - #endif - } - - /* Flip geometric normal with respect to the flip surface flag */ - flip_surface ^= ((struct geometry*)hit->prim.shape__)->flip_surface; - if(flip_surface) f3_minus(hit->normal, hit->normal); -} - static INLINE enum RTCBuildQuality accel_struct_quality_to_rtc_build_quality (enum s3d_accel_struct_quality quality) @@ -885,7 +783,7 @@ scene_view_compute_nprims_cdf nprims += 1; break; case GEOM_INSTANCE: - /* The instance CDF was computed during its scnview synchronisation */ + /* The instance CDF was computed during its scnview synchronisation */ len = darray_nprims_cdf_size_get (&geom->data.instance->scnview->nprims_cdf); if(len) { @@ -1251,120 +1149,6 @@ s3d_scene_view_get_mask(struct s3d_scene_view* scnview, int* mask) } res_T -s3d_scene_view_trace_ray - (struct s3d_scene_view* scnview, - const float org[3], - const float dir[3], - const float range[2], - void* ray_data, - struct s3d_hit* hit) -{ - struct RTCRayHit ray_hit; - struct intersect_context intersect_ctx; - size_t i; - - if(!scnview || !org || !dir || !range || !hit) - return RES_BAD_ARG; - if(!f3_is_normalized(dir)) { - log_error(scnview->scn->dev, - "%s: unnormalized ray direction {%g, %g, %g}.\n", - FUNC_NAME, SPLIT3(dir)); - return RES_BAD_ARG; - } - if(range[0] < 0) { - log_error(scnview->scn->dev, - "%s: invalid ray range [%g, %g] - it must be in [0, INF).\n", - FUNC_NAME, range[0], range[1]); - return RES_BAD_ARG; - } - if((scnview->mask & S3D_TRACE) == 0) { - log_error(scnview->scn->dev, - "%s: the S3D_TRACE flag is not active onto the submitted scene view.\n", - FUNC_NAME); - return RES_BAD_OP; - } - if(range[0] > range[1]) { /* Degenerated range <=> disabled ray */ - *hit = S3D_HIT_NULL; - return RES_OK; - } - - /* Initialise the ray */ - ray_hit.ray.org_x = org[0]; - ray_hit.ray.org_y = org[1]; - ray_hit.ray.org_z = org[2]; - ray_hit.ray.dir_x = dir[0]; - ray_hit.ray.dir_y = dir[1]; - ray_hit.ray.dir_z = dir[2]; - ray_hit.ray.tnear = range[0]; - ray_hit.ray.tfar = range[1]; - ray_hit.ray.time = FLT_MAX; /* Invalid fields */ - ray_hit.ray.mask = UINT_MAX; /* Invalid fields */ - ray_hit.ray.id = UINT_MAX; /* Invalid fields */ - ray_hit.ray.flags = UINT_MAX; /* Invalid fields */ - - /* Initialise the hit */ - ray_hit.hit.geomID = RTC_INVALID_GEOMETRY_ID; - FOR_EACH(i, 0, RTC_MAX_INSTANCE_LEVEL_COUNT) { - ray_hit.hit.instID[i] = RTC_INVALID_GEOMETRY_ID; - } - - /* Initialise the intersect context */ - rtcInitIntersectContext(&intersect_ctx.rtc); - intersect_ctx.ws_org[0] = org[0]; - intersect_ctx.ws_org[1] = org[1]; - intersect_ctx.ws_org[2] = org[2]; - intersect_ctx.ws_dir[0] = dir[0]; - intersect_ctx.ws_dir[1] = dir[1]; - intersect_ctx.ws_dir[2] = dir[2]; - intersect_ctx.scnview = scnview; - intersect_ctx.data = ray_data; - - /* Here we go! */ - rtcIntersect1(scnview->rtc_scn, &intersect_ctx.rtc, &ray_hit); - - hit_setup(scnview, &ray_hit, hit); - return RES_OK; -} - -res_T -s3d_scene_view_trace_rays - (struct s3d_scene_view* scnview, - const size_t nrays, - const int mask, - const float* origins, - const float* directions, - const float* ranges, - void* rays_data, - const size_t sizeof_ray_data, - struct s3d_hit* hits) -{ - size_t iray; - size_t iorg, idir, irange, idata; - size_t org_step, dir_step, range_step, data_step; - res_T res = RES_OK; - - if(!scnview) return RES_BAD_ARG; - if(!nrays) return RES_OK; - - org_step = mask & S3D_RAYS_SINGLE_ORIGIN ? 0 : 3; - dir_step = mask & S3D_RAYS_SINGLE_DIRECTION ? 0 : 3; - range_step = mask & S3D_RAYS_SINGLE_RANGE ? 0 : 2; - data_step = (mask & S3D_RAYS_SINGLE_DATA) || !rays_data ? 0 : sizeof_ray_data; - iorg = idir = irange = idata = 0; - - FOR_EACH(iray, 0, nrays) { - res = s3d_scene_view_trace_ray(scnview, origins+iorg, directions+idir, - ranges+irange, (char*)rays_data+idata, hits+iray); - if(UNLIKELY(res != RES_OK)) break; - iorg += org_step; - idir += dir_step; - irange += range_step; - idata += data_step; - } - return res; -} - -res_T s3d_scene_view_sample (struct s3d_scene_view* scnview, const float u, @@ -1770,37 +1554,3 @@ scene_view_destroy(struct s3d_scene_view* scnview) MEM_RM(scnview->scn->dev->allocator, scnview); } -/* Wrapper between an Embree and a Star-3D filter function */ -void -rtc_hit_filter_wrapper(const struct RTCFilterFunctionNArguments* args) -{ - struct s3d_hit hit; - struct RTCRayHit ray_hit; - struct intersect_context* ctx; - struct geometry* geom; - struct hit_filter* filter; - ASSERT(args && args->N == 1 && args->context && args->valid[0] != 0); - - rtc_rayN_get_ray(args->ray, args->N, 0, &ray_hit.ray); - rtc_hitN_get_hit(args->hit, args->N, 0, &ray_hit.hit); - - ctx = CONTAINER_OF(args->context, struct intersect_context, rtc); - - geom = args->geometryUserPtr; - switch(geom->type) { - case GEOM_MESH: - filter = &geom->data.mesh->filter; - break; - case GEOM_SPHERE: - filter = &geom->data.sphere->filter; - break; - default: FATAL("Unreachable code\n"); break; - } - ASSERT(filter->func); - - hit_setup(ctx->scnview, &ray_hit, &hit); - if(filter->func(&hit, ctx->ws_org, ctx->ws_dir, ctx->data, filter->data)) { - args->valid[0] = 0; - } -} - diff --git a/src/s3d_scene_view_closest_point.c b/src/s3d_scene_view_closest_point.c @@ -0,0 +1,440 @@ +/* Copyright (C) 2015-2019 |Meso|Star> (contact@meso-star.com) + * + * This software is a computer program whose purpose is to describe a + * virtual 3D environment that can be ray-traced and sampled both robustly + * and efficiently. + * + * This software is governed by the CeCILL license under French law and + * abiding by the rules of distribution of free software. You can use, + * modify and/or redistribute the software under the terms of the CeCILL + * license as circulated by CEA, CNRS and INRIA at the following URL + * "http://www.cecill.info". + * + * As a counterpart to the access to the source code and rights to copy, + * modify and redistribute granted by the license, users are provided only + * with a limited warranty and the software's author, the holder of the + * economic rights, and the successive licensors have only limited + * liability. + * + * In this respect, the user's attention is drawn to the risks associated + * with loading, using, modifying and/or developing or reproducing the + * software by the user in light of its specific status of free software, + * that may mean that it is complicated to manipulate, and that also + * therefore means that it is reserved for developers and experienced + * professionals having in-depth computer knowledge. Users are therefore + * encouraged to load and test the software's suitability as regards their + * requirements in conditions enabling the security of their systems and/or + * data to be ensured and, more generally, to use and operate it in the + * same conditions as regards security. + * + * The fact that you are presently reading this means that you have had + * knowledge of the CeCILL license and that you accept its terms. */ + +#include "s3d.h" +#include "s3d_device_c.h" +#include "s3d_instance.h" +#include "s3d_geometry.h" +#include "s3d_mesh.h" +#include "s3d_scene_view_c.h" +#include "s3d_sphere.h" + +#include <rsys/float3.h> +#include <rsys/float33.h> + +struct point_query_context { + struct RTCPointQueryContext rtc; + struct s3d_scene_view* scnview; + void* data; /* Per point query defined data */ +}; + +/******************************************************************************* + * Helper functions + ******************************************************************************/ +static INLINE float* +closest_point_triangle + (const float p[3], /* Point */ + const float a[3], /* 1st triangle vertex */ + const float b[3], /* 2nd triangle vertex */ + const float c[3], /* 3rd triangle vertex */ + float closest_pt[3], /* Closest position */ + float uv[2]) /* UV of the closest position */ +{ + float ab[3], ac[3], ap[3], bp[3], cp[3]; + float d1, d2, d3, d4, d5, d6; + float va, vb, vc; + float rcp_triangle_area; + float v, w; + ASSERT(p && a && b && c && closest_pt && uv); + + f3_sub(ab, b, a); + f3_sub(ac, c, a); + + /* Check if the closest point is the triangle vertex 'a' */ + f3_sub(ap, p, a); + d1 = f3_dot(ab, ap); + d2 = f3_dot(ac, ap); + if(d1 <= 0.f && d2 <= 0.f) { + uv[0] = 1.f; + uv[1] = 0.f; + return f3_set(closest_pt, a); + } + + /* Check if the closest point is the triangle vertex 'b' */ + f3_sub(bp, p, b); + d3 = f3_dot(ab, bp); + d4 = f3_dot(ac, bp); + if(d3 >= 0.f && d4 <= d3) { + uv[0] = 0.f; + uv[1] = 1.f; + return f3_set(closest_pt, b); + } + + /* Check if the closest point is the triangle vertex 'c' */ + f3_sub(cp, p, c); + d5 = f3_dot(ab, cp); + d6 = f3_dot(ac, cp); + if(d6 >= 0.f && d5 <= d6) { + uv[0] = 0.f; + uv[1] = 0.f; + return f3_set(closest_pt, c); + } + + /* Check if the closest point is on the triangle edge 'ab' */ + vc = d1*d4 - d3*d2; + if(vc <= 0.f && d1 >= 0.f && d3 <= 0.f) { + const float s = d1 / (d1 - d3); + closest_pt[0] = a[0] + s*ab[0]; + closest_pt[1] = a[1] + s*ab[1]; + closest_pt[2] = a[2] + s*ab[2]; + uv[0] = 1.f - s; + uv[1] = s; + return closest_pt; + } + + /* Check if the closest point is on the triangle edge 'ac' */ + vb = d5*d2 - d1*d6; + if(vb <= 0.f && d2 >= 0.f && d6 <= 0.f) { + const float s = d2 / (d2 - d6); + closest_pt[0] = a[0] + s*ac[0]; + closest_pt[1] = a[1] + s*ac[1]; + closest_pt[2] = a[2] + s*ac[2]; + uv[0] = 1.f - s; + uv[1] = 0.f; + return closest_pt; + } + + /* Check if the closest point is on the triangle edge 'bc' */ + va = d3*d6 - d5*d4; + if(va <= 0.f && (d4 - d3) >= 0.f && (d5 - d6) >= 0.f) { + const float s = (d4 - d3) / ((d4 - d3) + (d5 - d6)); + closest_pt[0] = b[0] + s*(c[0] - b[0]); + closest_pt[1] = b[1] + s*(c[1] - b[1]); + closest_pt[2] = b[2] + s*(c[2] - b[2]); + uv[0] = 0.f; + uv[1] = 1.f - s; + return closest_pt; + } + + /* The closest point lies in the triangle: compute its barycentric + * coordinates */ + rcp_triangle_area = 1.f / (va + vb + vc); + v = vb * rcp_triangle_area; + w = vc * rcp_triangle_area; + + /* Save the uv barycentric coordinates */ + uv[0] = 1.f - v - w; + uv[1] = v; + + ASSERT(eq_epsf(uv[0] + uv[1] + w, 1.f, 1.e-4f)); + + /* Use the barycentric coordinates to compute the world space position of the + * closest point onto the triangle */ + closest_pt[0] = a[0] + v*ab[0] + w*ac[0]; + closest_pt[1] = a[1] + v*ab[1] + w*ac[1]; + closest_pt[2] = a[2] + v*ab[2] + w*ac[2]; + return closest_pt; +} + +static bool +closest_point_mesh + (struct RTCPointQueryFunctionArguments* args, + struct geometry* geom, + struct geometry* inst, /* Can be NULL */ + void* query_data) +{ + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_hit* out_hit = NULL; + struct hit_filter* filter = NULL; + const uint32_t* ids = NULL; + float closest_point[3]; + float query_pos[3]; + float v0[3], v1[3], v2[3]; + float E0[3], E1[3], Ng[3]; + float vec[3]; + float uv[2]; + float dst; + int flip_surface = 0; + ASSERT(args && geom && geom->type == GEOM_MESH); + ASSERT(args->primID < mesh_get_ntris(geom->data.mesh)); + + /* Fetch triangle indices */ + ids = mesh_get_ids(geom->data.mesh) + args->primID*3/*#indices per triangle*/; + + /* Fetch triangle vertices */ + ASSERT(geom->data.mesh->attribs_type[S3D_POSITION] == S3D_FLOAT3); + f3_set(v0, mesh_get_pos(geom->data.mesh) + ids[0]*3/*#coords per vertex*/); + f3_set(v1, mesh_get_pos(geom->data.mesh) + ids[1]*3/*#coords per vertex*/); + f3_set(v2, mesh_get_pos(geom->data.mesh) + ids[2]*3/*#coords per vertex*/); + + /* Local copy of the query position */ + query_pos[0] = args->query->x; + query_pos[1] = args->query->y; + query_pos[2] = args->query->z; + + if(args->context->instStackSize) { /* The mesh is instantiated */ + const float* transform; + transform = inst->data.instance->transform; + ASSERT(args->context->instStackSize == 1); + ASSERT(args->similarityScale == 1); + ASSERT(inst && inst->type == GEOM_INSTANCE); + ASSERT(f3_eq_eps(transform+0, args->context->inst2world[0]+0, 1.e-6f)); + ASSERT(f3_eq_eps(transform+3, args->context->inst2world[0]+4, 1.e-6f)); + ASSERT(f3_eq_eps(transform+6, args->context->inst2world[0]+8, 1.e-6f)); + ASSERT(f3_eq_eps(transform+9, args->context->inst2world[0]+12,1.e-6f)); + + /* Transform the triangle vertices in world space */ + f3_add(v0, f33_mulf3(v0, transform, v0), transform+9); + f3_add(v1, f33_mulf3(v1, transform, v1), transform+9); + f3_add(v2, f33_mulf3(v2, transform, v2), transform+9); + + flip_surface = inst->flip_surface; + } + + /* Compute the closest point onto the triangle from the submitted point */ + closest_point_triangle(query_pos, v0, v1, v2, closest_point, uv); + + f3_sub(vec, closest_point, query_pos); + dst = f3_len(vec); + if(dst >= args->query->radius) return 0; + + /* Compute the triangle normal in world space */ + f3_sub(E0, v2, v0); + f3_sub(E1, v1, v0); + f3_cross(Ng, E0, E1); + + /* Flip the geometric normal wrt the flip surface flag */ + flip_surface ^= geom->flip_surface; + if(flip_surface) f3_minus(Ng, Ng); + + /* Setup the hit */ + hit.prim.shape__ = geom; + hit.prim.inst__ = inst; + hit.distance = dst; + hit.uv[0] = uv[0]; + hit.uv[1] = uv[1]; + hit.normal[0] = Ng[0]; + hit.normal[1] = Ng[1]; + hit.normal[2] = Ng[2]; + hit.prim.prim_id = args->primID; + hit.prim.geom_id = geom->name; + hit.prim.inst_id = inst ? inst->name : S3D_INVALID_ID; + hit.prim.scene_prim_id = + hit.prim.prim_id + + geom->scene_prim_id_offset + + (inst ? inst->scene_prim_id_offset : 0); + + /* `vec' is the direction along which the closest point was found. We thus + * submit it to the filter function as the direction corresponding to the + * computed hit */ + filter = &geom->data.mesh->filter; + if(filter->func + && filter->func(&hit, query_pos, vec, query_data, filter->data)) { + return 0; /* This point is filtered. Discard it! */ + } + + /* Update output data */ + out_hit = args->userPtr; + *out_hit = hit; + + /* Shrink the query radius */ + args->query->radius = dst; + + return 1; /* Notify that the query radius was updated */ +} + +static bool +closest_point_sphere + (struct RTCPointQueryFunctionArguments* args, + struct geometry* geom, + struct geometry* inst, + void* query_data) +{ + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_hit* out_hit = NULL; + struct hit_filter* filter = NULL; + float query_pos[3]; + float sphere_pos[3]; + float Ng[3]; + float dir[3]; + float uv[2]; + float dst; + float len; + int flip_surface = 0; + ASSERT(args && geom && geom->type == GEOM_SPHERE && args->primID == 0); + + /* Local copy of the query position */ + query_pos[0] = args->query->x; + query_pos[1] = args->query->y; + query_pos[2] = args->query->z; + + f3_set(sphere_pos, geom->data.sphere->pos); + if(args->context->instStackSize) { /* The sphere is instantiated */ + const float* transform; + transform = inst->data.instance->transform; + ASSERT(args->context->instStackSize == 1); + ASSERT(args->similarityScale == 1); + ASSERT(inst && inst->type == GEOM_INSTANCE); + ASSERT(f3_eq_eps(transform+0, args->context->inst2world[0]+0, 1.e-6f)); + ASSERT(f3_eq_eps(transform+3, args->context->inst2world[0]+4, 1.e-6f)); + ASSERT(f3_eq_eps(transform+6, args->context->inst2world[0]+8, 1.e-6f)); + ASSERT(f3_eq_eps(transform+9, args->context->inst2world[0]+12,1.e-6f)); + + /* Transform the sphere position in world space */ + f33_mulf3(sphere_pos, transform, sphere_pos); + f3_add(sphere_pos, transform+9, sphere_pos); + + flip_surface = inst->flip_surface; + } + + /* Compute the distance from the the query pos to the sphere center */ + f3_sub(Ng, query_pos, sphere_pos); + len = f3_len(Ng); + + /* Evaluate the distance from the query pos to the sphere surface */ + dst = fabsf(len - geom->data.sphere->radius); + + /* The closest point onto the sphere is outside the query radius */ + if(dst >= args->query->radius) + return 0; + + /* Compute the uv of the found point */ + f3_divf(Ng, Ng, len); /* Normalize the hit normal */ + sphere_normal_to_uv(Ng, uv); + + /* Flip the geometric normal wrt the flip surface flag */ + flip_surface ^= geom->flip_surface; + if(flip_surface) f3_minus(Ng, Ng); + + /* Setup the hit */ + hit.prim.shape__ = geom; + hit.prim.inst__ = inst; + hit.distance = dst; + hit.uv[0] = uv[0]; + hit.uv[1] = uv[1]; + hit.normal[0] = Ng[0]; + hit.normal[1] = Ng[1]; + hit.normal[2] = Ng[2]; + hit.prim.prim_id = args->primID; + hit.prim.geom_id = geom->name; + hit.prim.inst_id = inst ? inst->name : S3D_INVALID_ID; + hit.prim.scene_prim_id = + hit.prim.prim_id + + geom->scene_prim_id_offset + + (inst ? inst->scene_prim_id_offset : 0); + + /* Use the reversed geometric normal as the hit direction since it is along + * this vector that the closest point was effectively computed */ + f3_minus(dir, Ng); + filter = &geom->data.sphere->filter; + if(filter->func + && filter->func(&hit, query_pos, dir, query_data, filter->data)) { + return 0; + } + + /* Update output data */ + out_hit = args->userPtr; + *out_hit = hit; + + /* Shrink the query radius */ + args->query->radius = dst; + + return 1; /* Notify that the query radius was updated */ +} + +static bool +closest_point(struct RTCPointQueryFunctionArguments* args) +{ + struct point_query_context* ctx = NULL; + struct geometry* geom = NULL; + struct geometry* inst = NULL; + bool query_radius_is_upd = false; + ASSERT(args); + + ctx = CONTAINER_OF(args->context, struct point_query_context, rtc); + if(args->context->instStackSize == 0) { + geom = scene_view_geometry_from_embree_id + (ctx->scnview, args->geomID); + } else { + ASSERT(args->context->instStackSize == 1); + ASSERT(args->similarityScale == 1); + inst = scene_view_geometry_from_embree_id + (ctx->scnview, args->context->instID[0]); + geom = scene_view_geometry_from_embree_id + (inst->data.instance->scnview, args->geomID); + } + + switch(geom->type) { + case GEOM_MESH: + query_radius_is_upd = closest_point_mesh(args, geom, inst, ctx->data); + break; + case GEOM_SPHERE: + query_radius_is_upd = closest_point_sphere(args, geom, inst, ctx->data); + break; + default: FATAL("Unreachable code\n"); break; + } + return query_radius_is_upd; +} + +/******************************************************************************* + * Exported functions + ******************************************************************************/ +res_T +s3d_scene_view_closest_point + (struct s3d_scene_view* scnview, + const float pos[3], + const float radius, + void* query_data, + struct s3d_hit* hit) +{ + struct RTCPointQuery query; + struct point_query_context query_ctx; + + if(!scnview || !pos || radius <= 0 || !hit) + return RES_BAD_ARG; + if((scnview->mask & S3D_TRACE) == 0) { + log_error(scnview->scn->dev, + "%s: the S3D_TRACE flag is not active onto the submitted scene view.\n", + FUNC_NAME); + return RES_BAD_OP; + } + + *hit = S3D_HIT_NULL; + + /* Initialise the point query */ + query.x = pos[0]; + query.y = pos[1]; + query.z = pos[2]; + query.radius = radius; + query.time = FLT_MAX; /* Invalid fields */ + + /* Initialise the point query context */ + rtcInitPointQueryContext(&query_ctx.rtc); + query_ctx.scnview = scnview; + query_ctx.data = query_data; + + /* Here we go! */ + rtcPointQuery(scnview->rtc_scn, &query, &query_ctx.rtc, closest_point, hit); + + return RES_OK; +} + diff --git a/src/s3d_scene_view_trace_ray.c b/src/s3d_scene_view_trace_ray.c @@ -0,0 +1,302 @@ +/* Copyright (C) 2015-2019 |Meso|Star> (contact@meso-star.com) + * + * This software is a computer program whose purpose is to describe a + * virtual 3D environment that can be ray-traced and sampled both robustly + * and efficiently. + * + * This software is governed by the CeCILL license under French law and + * abiding by the rules of distribution of free software. You can use, + * modify and/or redistribute the software under the terms of the CeCILL + * license as circulated by CEA, CNRS and INRIA at the following URL + * "http://www.cecill.info". + * + * As a counterpart to the access to the source code and rights to copy, + * modify and redistribute granted by the license, users are provided only + * with a limited warranty and the software's author, the holder of the + * economic rights, and the successive licensors have only limited + * liability. + * + * In this respect, the user's attention is drawn to the risks associated + * with loading, using, modifying and/or developing or reproducing the + * software by the user in light of its specific status of free software, + * that may mean that it is complicated to manipulate, and that also + * therefore means that it is reserved for developers and experienced + * professionals having in-depth computer knowledge. Users are therefore + * encouraged to load and test the software's suitability as regards their + * requirements in conditions enabling the security of their systems and/or + * data to be ensured and, more generally, to use and operate it in the + * same conditions as regards security. + * + * The fact that you are presently reading this means that you have had + * knowledge of the CeCILL license and that you accept its terms. */ + +#include "s3d.h" +#include "s3d_c.h" +#include "s3d_device_c.h" +#include "s3d_instance.h" +#include "s3d_geometry.h" +#include "s3d_mesh.h" +#include "s3d_sphere.h" +#include "s3d_scene_view_c.h" + +#include <rsys/float33.h> +#include <limits.h> + +struct intersect_context { + struct RTCIntersectContext rtc; + struct s3d_scene_view* scnview; + void* data; /* Per ray user defined data */ + float ws_org[3]; /* World space ray origin */ + float ws_dir[3]; /* World space ray direction */ +}; + +/******************************************************************************* + * Helper functions + ******************************************************************************/ +static INLINE void +hit_setup + (struct s3d_scene_view* scnview, + const struct RTCRayHit* ray_hit, + struct s3d_hit* hit) +{ + float w; + char flip_surface = 0; + + ASSERT(scnview && hit && ray_hit); + + if(ray_hit->hit.geomID == RTC_INVALID_GEOMETRY_ID) { /* No hit */ + *hit = S3D_HIT_NULL; + return; + } + + hit->normal[0] = ray_hit->hit.Ng_x; + hit->normal[1] = ray_hit->hit.Ng_y; + hit->normal[2] = ray_hit->hit.Ng_z; + hit->distance = ray_hit->ray.tfar; + + if(ray_hit->hit.instID[0] == RTC_INVALID_GEOMETRY_ID) { + struct geometry* geom_shape; + geom_shape = scene_view_geometry_from_embree_id(scnview, ray_hit->hit.geomID); + hit->prim.shape__ = geom_shape; + hit->prim.inst__ = NULL; + hit->prim.prim_id = ray_hit->hit.primID; + hit->prim.geom_id = geom_shape->name; + hit->prim.inst_id = S3D_INVALID_ID; + hit->prim.scene_prim_id = /* Compute the "scene space" primitive id */ + hit->prim.prim_id /* Mesh space */ + + geom_shape->scene_prim_id_offset; /* Scene space */ + + } else { /* The hit shape is instantiated */ + /* Retrieve the hit instance */ + struct geometry* geom_inst; + struct geometry* geom_shape; + float transform[9]; + geom_inst = scene_view_geometry_from_embree_id + (scnview, ray_hit->hit.instID[0]); + geom_shape = scene_view_geometry_from_embree_id + (geom_inst->data.instance->scnview, ray_hit->hit.geomID); + hit->prim.shape__ = geom_shape; + hit->prim.inst__ = geom_inst; + hit->prim.prim_id = ray_hit->hit.primID; + hit->prim.geom_id = geom_shape->name; + hit->prim.inst_id = geom_inst->name; + hit->prim.scene_prim_id = /* Compute the "scene space" primitive id */ + hit->prim.prim_id /* Shape space */ + + geom_shape->scene_prim_id_offset /* Inst space */ + + geom_inst->scene_prim_id_offset; /* Scene space */ + + flip_surface = geom_inst->flip_surface; + ASSERT(hit->prim.inst__); + ASSERT(((struct geometry*)hit->prim.inst__)->type == GEOM_INSTANCE); + + /* Transform the normal in world space */ + f33_invtrans(transform, geom_inst->data.instance->transform); + f33_mulf3(hit->normal, transform, hit->normal); + } + ASSERT(hit->prim.shape__); + ASSERT(((struct geometry*)hit->prim.shape__)->type == GEOM_MESH + ||((struct geometry*)hit->prim.shape__)->type == GEOM_SPHERE); + + hit->uv[0] = ray_hit->hit.u; + hit->uv[1] = ray_hit->hit.v; + + if(((struct geometry*)hit->prim.shape__)->type == GEOM_MESH) { + w = 1.f - hit->uv[0] - hit->uv[1]; + ASSERT(w <= 1.f); /* This may not occurs */ + if(w < 0.f) { /* Handle precision error */ + if(hit->uv[0] > hit->uv[1]) hit->uv[0] += w; + else hit->uv[1] += w; + w = 0.f; + } + + /* Embree stores on the u and v ray parameters the barycentric coordinates of + * the hit with respect to the second and third triangle vertices, + * respectively. The following code computes the barycentric coordinates of + * the hit for the first and second triangle vertices */ + hit->uv[1] = hit->uv[0]; + hit->uv[0] = w; + + /* In Embree3 the normal orientation is flipped wrt to Star-3D convention */ + #if RTC_VERSION_MAJOR >= 3 + f3_minus(hit->normal, hit->normal); + #endif + } + + /* Flip geometric normal with respect to the flip surface flag */ + flip_surface ^= ((struct geometry*)hit->prim.shape__)->flip_surface; + if(flip_surface) f3_minus(hit->normal, hit->normal); +} + +/******************************************************************************* + * Exported functions + ******************************************************************************/ +res_T +s3d_scene_view_trace_ray + (struct s3d_scene_view* scnview, + const float org[3], + const float dir[3], + const float range[2], + void* ray_data, + struct s3d_hit* hit) +{ + struct RTCRayHit ray_hit; + struct intersect_context intersect_ctx; + size_t i; + + if(!scnview || !org || !dir || !range || !hit) + return RES_BAD_ARG; + if(!f3_is_normalized(dir)) { + log_error(scnview->scn->dev, + "%s: unnormalized ray direction {%g, %g, %g}.\n", + FUNC_NAME, SPLIT3(dir)); + return RES_BAD_ARG; + } + if(range[0] < 0) { + log_error(scnview->scn->dev, + "%s: invalid ray range [%g, %g] - it must be in [0, INF).\n", + FUNC_NAME, range[0], range[1]); + return RES_BAD_ARG; + } + if((scnview->mask & S3D_TRACE) == 0) { + log_error(scnview->scn->dev, + "%s: the S3D_TRACE flag is not active onto the submitted scene view.\n", + FUNC_NAME); + return RES_BAD_OP; + } + if(range[0] > range[1]) { /* Degenerated range <=> disabled ray */ + *hit = S3D_HIT_NULL; + return RES_OK; + } + + /* Initialise the ray */ + ray_hit.ray.org_x = org[0]; + ray_hit.ray.org_y = org[1]; + ray_hit.ray.org_z = org[2]; + ray_hit.ray.dir_x = dir[0]; + ray_hit.ray.dir_y = dir[1]; + ray_hit.ray.dir_z = dir[2]; + ray_hit.ray.tnear = range[0]; + ray_hit.ray.tfar = range[1]; + ray_hit.ray.time = FLT_MAX; /* Invalid fields */ + ray_hit.ray.mask = UINT_MAX; /* Invalid fields */ + ray_hit.ray.id = UINT_MAX; /* Invalid fields */ + ray_hit.ray.flags = UINT_MAX; /* Invalid fields */ + + /* Initialise the hit */ + ray_hit.hit.geomID = RTC_INVALID_GEOMETRY_ID; + FOR_EACH(i, 0, RTC_MAX_INSTANCE_LEVEL_COUNT) { + ray_hit.hit.instID[i] = RTC_INVALID_GEOMETRY_ID; + } + + /* Initialise the intersect context */ + rtcInitIntersectContext(&intersect_ctx.rtc); + intersect_ctx.ws_org[0] = org[0]; + intersect_ctx.ws_org[1] = org[1]; + intersect_ctx.ws_org[2] = org[2]; + intersect_ctx.ws_dir[0] = dir[0]; + intersect_ctx.ws_dir[1] = dir[1]; + intersect_ctx.ws_dir[2] = dir[2]; + intersect_ctx.scnview = scnview; + intersect_ctx.data = ray_data; + + /* Here we go! */ + rtcIntersect1(scnview->rtc_scn, &intersect_ctx.rtc, &ray_hit); + + hit_setup(scnview, &ray_hit, hit); + return RES_OK; +} + +res_T +s3d_scene_view_trace_rays + (struct s3d_scene_view* scnview, + const size_t nrays, + const int mask, + const float* origins, + const float* directions, + const float* ranges, + void* rays_data, + const size_t sizeof_ray_data, + struct s3d_hit* hits) +{ + size_t iray; + size_t iorg, idir, irange, idata; + size_t org_step, dir_step, range_step, data_step; + res_T res = RES_OK; + + if(!scnview) return RES_BAD_ARG; + if(!nrays) return RES_OK; + + org_step = mask & S3D_RAYS_SINGLE_ORIGIN ? 0 : 3; + dir_step = mask & S3D_RAYS_SINGLE_DIRECTION ? 0 : 3; + range_step = mask & S3D_RAYS_SINGLE_RANGE ? 0 : 2; + data_step = (mask & S3D_RAYS_SINGLE_DATA) || !rays_data ? 0 : sizeof_ray_data; + iorg = idir = irange = idata = 0; + + FOR_EACH(iray, 0, nrays) { + res = s3d_scene_view_trace_ray(scnview, origins+iorg, directions+idir, + ranges+irange, (char*)rays_data+idata, hits+iray); + if(UNLIKELY(res != RES_OK)) break; + iorg += org_step; + idir += dir_step; + irange += range_step; + idata += data_step; + } + return res; +} + +/******************************************************************************* + * Local functions + ******************************************************************************/ +/* Wrapper between an Embree and a Star-3D filter function */ +void +rtc_hit_filter_wrapper(const struct RTCFilterFunctionNArguments* args) +{ + struct s3d_hit hit; + struct RTCRayHit ray_hit; + struct intersect_context* ctx; + struct geometry* geom; + struct hit_filter* filter; + ASSERT(args && args->N == 1 && args->context && args->valid[0] != 0); + + rtc_rayN_get_ray(args->ray, args->N, 0, &ray_hit.ray); + rtc_hitN_get_hit(args->hit, args->N, 0, &ray_hit.hit); + + ctx = CONTAINER_OF(args->context, struct intersect_context, rtc); + + geom = args->geometryUserPtr; + switch(geom->type) { + case GEOM_MESH: + filter = &geom->data.mesh->filter; + break; + case GEOM_SPHERE: + filter = &geom->data.sphere->filter; + break; + default: FATAL("Unreachable code\n"); break; + } + ASSERT(filter->func); + + hit_setup(ctx->scnview, &ray_hit, &hit); + if(filter->func(&hit, ctx->ws_org, ctx->ws_dir, ctx->data, filter->data)) { + args->valid[0] = 0; + } +} diff --git a/src/s3d_sphere.h b/src/s3d_sphere.h @@ -100,4 +100,26 @@ sphere_compute_volume(const struct sphere* sphere) return (float)(4*PI*r*r*r/3); } +static FINLINE void +sphere_normal_to_uv(const float normal[3], float uv[2]) +{ + float u, v, cos_theta; + ASSERT(normal && uv && f3_is_normalized(normal)); + + cos_theta = normal[2]; + + v = (1.f - cos_theta) * 0.5f; + if(absf(cos_theta) == 1) { + u = 0; + } else if(eq_epsf(normal[0], 0.f, 1.e-6f)) { + u = normal[1] > 0 ? 0.25f : 0.75f; + } else { + double phi = atan2f(normal[1], normal[0]); /* phi in [-PI, PI] */ + if(phi < 0) phi = 2*PI + phi; /* phi in [0, 2PI] */ + u = (float)(phi / (2*PI)); + } + uv[0] = u; + uv[1] = v; +} + #endif /* S3D_SPHERE_H */ diff --git a/src/test_s3d_closest_point.c b/src/test_s3d_closest_point.c @@ -0,0 +1,822 @@ +/* Copyright (C) 2015-2019 |Meso|Star> (contact@meso-star.com) + * + * This software is a computer program whose purpose is to describe a + * virtual 3D environment that can be ray-traced and sampled both robustly + * and efficiently. + * + * This software is governed by the CeCILL license under French law and + * abiding by the rules of distribution of free software. You can use, + * modify and/or redistribute the software under the terms of the CeCILL + * license as circulated by CEA, CNRS and INRIA at the following URL + * "http://www.cecill.info". + * + * As a counterpart to the access to the source code and rights to copy, + * modify and redistribute granted by the license, users are provided only + * with a limited warranty and the software's author, the holder of the + * economic rights, and the successive licensors have only limited + * liability. + * + * In this respect, the user's attention is drawn to the risks associated + * with loading, using, modifying and/or developing or reproducing the + * software by the user in light of its specific status of free software, + * that may mean that it is complicated to manipulate, and that also + * therefore means that it is reserved for developers and experienced + * professionals having in-depth computer knowledge. Users are therefore + * encouraged to load and test the software's suitability as regards their + * requirements in conditions enabling the security of their systems and/or + * data to be ensured and, more generally, to use and operate it in the + * same conditions as regards security. + * + * The fact that you are presently reading this means that you have had + * knowledge of the CeCILL license and that you accept its terms. */ + +#include "s3d.h" +#include "test_s3d_cbox.h" +#include "test_s3d_utils.h" + +#include <rsys/float3.h> +#include <limits.h> + +#define ON_EDGE_EPSILON 1.e-4f +#define POSITION_EPSILON 1.e-3f + +struct closest_pt { + float pos[3]; + float normal[3]; + float dst; + unsigned iprim; + unsigned igeom; + unsigned iinst; +}; + +#define CLOSEST_PT_NULL__ { \ + {0,0,0}, \ + {0,0,0}, \ + FLT_MAX, \ + S3D_INVALID_ID, \ + S3D_INVALID_ID, \ + S3D_INVALID_ID \ +} + +static const struct closest_pt CLOSEST_PT_NULL = CLOSEST_PT_NULL__; + +/******************************************************************************* + * Helper functions + ******************************************************************************/ +/* Function that computes the point onto the triangle that ensures the minimum + * distance between the submitted `pos' and the triangle. Use this routine to + * cross check the result of the s3d_scene_view_closest_point function that + * internally relies on a more efficient implementation */ +static float* +closest_point_triangle + (const float p[3], /* Input pos */ + const float a[3], /* 1st triangle vertex */ + const float b[3], /* 2nd triangle vertex */ + const float c[3], /* 3rd triangle vertex */ + float pt[3]) /* Closest point of pos onto the triangle */ +{ + float N[3]; /* Triangle normal */ + float Nab[3], Nbc[3], Nca[3]; /* Edge normals */ + float ab[3], ac[3], bc[3]; + float ap[3], bp[3], cp[3]; + float d1, d2, d3, d4, d5, d6, d; + CHK(p && a && b && c && pt); + + f3_normalize(ab, f3_sub(ab, b, a)); + f3_normalize(ac, f3_sub(ac, c, a)); + f3_normalize(bc, f3_sub(bc, c, b)); + + /* Compute the triangle normal */ + f3_cross(N, ac, ab); + + /* Check if the nearest point is the 1st triangle vertex */ + f3_sub(ap, p, a); + d1 = f3_dot(ab, ap); + d2 = f3_dot(ac, ap); + if(d1 <= 0 && d2 <= 0) return f3_set(pt, a); + + /* Check if the nearest point is the 2nd triangle vertex */ + f3_sub(bp, p, b); + d3 = f3_dot(bc, bp); + d4 =-f3_dot(ab, bp); + if(d3 <= 0 && d4 <= 0) return f3_set(pt, b); + + /* Check if the nearest point is the 3rd triangle vertex */ + f3_sub(cp, p, c); + d5 =-f3_dot(ac, cp); + d6 =-f3_dot(bc, cp); + if(d5 <= 0 && d6 <= 0) return f3_set(pt, c); + + /* Check if the nearest point is on the 1st triangle edge */ + f3_normalize(Nbc, f3_cross(Nab, ab, N)); + if(f3_dot(Nab, ap) <= 0) { + return f3_add(pt, a, f3_mulf(pt, ab, d1)); + } + + /* Check if the nearest point is on the 2nd triangle edge */ + f3_normalize(Nbc, f3_cross(Nbc, bc, N)); + if(f3_dot(Nbc, bp) <= 0) { + return f3_add(pt, b, f3_mulf(pt, bc, d3)); + } + + /* Check if the nearest point is on the 3rd triangle edge */ + f3_normalize(Nca, f3_cross(Nca, ac, N)); + f3_minus(Nca, Nca); + if(f3_dot(Nca, cp) <= 0) { + return f3_add(pt, c, f3_mulf(pt, ac,-d5)); + } + + /* The nearest point is in the triangle */ + f3_normalize(N, N); + d = f3_dot(N, ap); + return f3_add(pt, p, f3_mulf(pt, N, -d)); +} + +static void +closest_point_mesh + (const float pos[3], + const float* verts, + const unsigned* ids, + const unsigned ntris, + const unsigned geom_id, + const unsigned inst_id, + struct closest_pt* pt) +{ + unsigned itri; + CHK(pos && verts && ids && pt); + + *pt = CLOSEST_PT_NULL; + pt->igeom = geom_id; + pt->iinst = inst_id; + pt->dst = FLT_MAX; + + /* Find the closest point on the mesh */ + FOR_EACH(itri, 0, ntris) { + float v0[3]; + float v1[3]; + float v2[3]; + float closest_pt[3]; + float vec[3]; + float dst; + + f3_set(v0, verts+ids[itri*3+0]*3); + f3_set(v1, verts+ids[itri*3+1]*3); + f3_set(v2, verts+ids[itri*3+2]*3); + + closest_point_triangle(pos, v0, v1, v2, closest_pt); + dst = f3_len(f3_sub(vec, closest_pt, pos)); + + if(dst < pt->dst) { + float E0[3], E1[3]; + f3_set(pt->pos, closest_pt); + pt->dst = dst; + pt->iprim = itri; + f3_sub(E0, v1, v0); + f3_sub(E1, v2, v0); + f3_cross(pt->normal, E1, E0); + f3_normalize(pt->normal, pt->normal); + } + } +} + +static void +closest_point_sphere + (const float pos[3], + const float sphere_org[3], + const float sphere_radius, + const unsigned geom_id, + const unsigned inst_id, + struct closest_pt* pt) +{ + float vec[3]; + float len; + CHK(pos && sphere_org && sphere_radius > 0 && pt); + + f3_sub(vec, pos, sphere_org); + len = f3_normalize(vec, vec); + + pt->dst = (float)fabs(len - sphere_radius); + f3_set(pt->normal, vec); + f3_add(pt->pos, sphere_org, f3_mulf(pt->pos, vec, sphere_radius)); + pt->iprim = 0; + pt->igeom = geom_id; + pt->iinst = inst_id; +} + +/* Check that `hit' roughly lies on an edge. */ +static int +hit_on_edge(const struct s3d_hit* hit) +{ + struct s3d_attrib v0, v1, v2, pos; + float E0[3], E1[3], N[3]; + float tri_2area; + float hit_2area0; + float hit_2area1; + float hit_2area2; + float hit_pos[3]; + + CHK(hit && !S3D_HIT_NONE(hit)); + + /* Retrieve the triangle vertices */ + CHK(s3d_triangle_get_vertex_attrib(&hit->prim, 0, S3D_POSITION, &v0)==RES_OK); + CHK(s3d_triangle_get_vertex_attrib(&hit->prim, 1, S3D_POSITION, &v1)==RES_OK); + CHK(s3d_triangle_get_vertex_attrib(&hit->prim, 2, S3D_POSITION, &v2)==RES_OK); + + /* Compute the triangle area * 2 */ + f3_sub(E0, v1.value, v0.value); + f3_sub(E1, v2.value, v0.value); + tri_2area = f3_len(f3_cross(N, E0, E1)); + + /* Compute the hit position */ + CHK(s3d_primitive_get_attrib(&hit->prim, S3D_POSITION, hit->uv, &pos) == RES_OK); + f3_set(hit_pos, pos.value); + + /* Compute areas */ + f3_sub(E0, v0.value, hit_pos); + f3_sub(E1, v1.value, hit_pos); + hit_2area0 = f3_len(f3_cross(N, E0, E1)); + f3_sub(E0, v1.value, hit_pos); + f3_sub(E1, v2.value, hit_pos); + hit_2area1 = f3_len(f3_cross(N, E0, E1)); + f3_sub(E0, v2.value, hit_pos); + f3_sub(E1, v0.value, hit_pos); + hit_2area2 = f3_len(f3_cross(N, E0, E1)); + + if(hit_2area0 / tri_2area < ON_EDGE_EPSILON + || hit_2area1 / tri_2area < ON_EDGE_EPSILON + || hit_2area2 / tri_2area < ON_EDGE_EPSILON) + return 1; + + return 0; +} + +static void +check_closest_point + (const struct s3d_hit* hit, + const struct closest_pt* pt, + const int hit_triangle) /* Define if `hit' lies on a triangle */ +{ + struct s3d_attrib attr; + float N[3]; + + CHK(hit && pt); + CHK(!S3D_HIT_NONE(hit)); + + CHK(s3d_primitive_get_attrib + (&hit->prim, S3D_POSITION, hit->uv, &attr) == RES_OK); + f3_normalize(N, hit->normal); + + if(!hit_triangle || !hit_on_edge(hit)) { + CHK(hit->prim.prim_id == pt->iprim); + } + + if(hit->prim.prim_id == pt->iprim + && hit->prim.geom_id == pt->igeom + && hit->prim.inst_id == pt->iinst) { + /* Due to numerical inaccuracies and/or the arbitrary order in which + * primitives are treated, 2 points on different primitive can have the + * same distance from the query position while their respective + * coordinates are not equal wrt POSITION_EPSILON. To avoid wrong + * assertion, we thus check the position returned by Star-3D against the + * manually computed position only if these positions lies on the same + * primitive */ + CHK(f3_eq_eps(pt->pos, attr.value, POSITION_EPSILON)); + CHK(f3_eq_eps(pt->normal, N, 1.e-4f)); + } + CHK(eq_epsf(hit->distance, pt->dst, 1.e-3f)); +} + +/******************************************************************************* + * Cornell box and sphere test + ******************************************************************************/ +struct instance { + float translation[3]; + unsigned id; +}; + +static void +check_closest_point_cbox_sphere + (const float pos[3], + const float sphere_org[3], + const float sphere_radius, + const unsigned walls_id, + const unsigned sphere_id, + const struct instance* instances, + const size_t ninstances, + struct s3d_hit* hit) +{ + struct closest_pt pt_walls = CLOSEST_PT_NULL; + struct closest_pt pt_sphere = CLOSEST_PT_NULL; + const struct closest_pt* pt = NULL; + CHK(pos && hit); + + if(!ninstances) { + closest_point_mesh(pos, cbox_walls, cbox_walls_ids, cbox_walls_ntris, + walls_id, S3D_INVALID_ID, &pt_walls); + closest_point_sphere(pos, sphere_org, sphere_radius, sphere_id, + S3D_INVALID_ID, &pt_sphere); + } else { + size_t iinst; + + pt_walls.dst = FLT_MAX; + FOR_EACH(iinst, 0, ninstances) { + struct closest_pt pt_walls_tmp; + struct closest_pt pt_sphere_tmp; + float pos_instance_space[3]; + + /* Transform query position in instance space */ + f3_sub(pos_instance_space, pos, instances[iinst].translation); + + closest_point_mesh(pos_instance_space, cbox_walls, cbox_walls_ids, + cbox_walls_ntris, walls_id, instances[iinst].id, &pt_walls_tmp); + closest_point_sphere(pos_instance_space, sphere_org, sphere_radius, + sphere_id, instances[iinst].id, &pt_sphere_tmp); + + if(pt_walls_tmp.dst < pt_walls.dst) { + pt_walls = pt_walls_tmp; + /* Transform query closest point in world space */ + f3_add(pt_walls.pos, pt_walls.pos, instances[iinst].translation); + } + if(pt_sphere_tmp.dst < pt_sphere.dst) { + pt_sphere = pt_sphere_tmp; + /* Transform query closest point in world space */ + f3_add(pt_sphere.pos, pt_sphere.pos, instances[iinst].translation); + } + } + } + + if(pt_walls.dst< pt_sphere.dst) { + pt = &pt_walls; + } else { + pt = &pt_sphere; + } + + check_closest_point(hit, pt, hit->prim.geom_id == walls_id); +} + +static void +test_cbox_sphere(struct s3d_device* dev) +{ + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL; + struct s3d_scene* scn = NULL; + struct s3d_shape* walls = NULL; + struct s3d_shape* sphere = NULL; + struct s3d_shape* inst0 = NULL; + struct s3d_shape* inst1 = NULL; + struct s3d_scene_view* scnview = NULL; + struct instance instances[2]; + struct cbox_desc cbox_desc; + size_t i; + float low[3], upp[3], mid[3], sz[3]; + float pos[3]; + float sphere_org[3]; + float sphere_radius; + unsigned walls_id, sphere_id; + + CHK(s3d_scene_create(dev, &scn) == RES_OK); + + /* Setup the cornell box walls */ + vdata.usage = S3D_POSITION; + vdata.type = S3D_FLOAT3; + vdata.get = cbox_get_position; + cbox_desc.vertices = cbox_walls; + cbox_desc.indices = cbox_walls_ids; + CHK(s3d_shape_create_mesh(dev, &walls) == RES_OK); + CHK(s3d_shape_get_id(walls, &walls_id) == RES_OK); + CHK(s3d_scene_attach_shape(scn, walls) == RES_OK); + CHK(s3d_mesh_setup_indexed_vertices(walls, cbox_walls_ntris, cbox_get_ids, + cbox_walls_nverts, &vdata, 1, &cbox_desc) == RES_OK); + + /* Compute the Cornell box AABB */ + CHK(s3d_scene_view_create(scn, S3D_GET_PRIMITIVE, &scnview) == RES_OK); + CHK(s3d_scene_view_get_aabb(scnview, low, upp) == RES_OK); + CHK(s3d_scene_view_ref_put(scnview) == RES_OK); + + /* Setup the sphere at the center of the cornell box */ + f3_mulf(mid, f3_add(mid, low, upp), 0.5f); + f3_sub(sz, upp, low); + f3_set(sphere_org, mid); + sphere_radius = MMIN(MMIN(sz[0], sz[1]), sz[2]) * 0.125f; /* 1/8 of the box */ + CHK(s3d_shape_create_sphere(dev, &sphere) == RES_OK); + CHK(s3d_shape_get_id(sphere, &sphere_id) == RES_OK); + CHK(s3d_scene_attach_shape(scn, sphere) == RES_OK); + CHK(s3d_sphere_setup(sphere, sphere_org, sphere_radius) == RES_OK); + + CHK(s3d_scene_view_create(scn, S3D_TRACE, &scnview) == RES_OK); + + /* Check point query on the scene */ + FOR_EACH(i, 0, 10000) { + /* Randomly generate a point in a bounding box that is 2 times the size of + * the scene AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]); + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]); + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]); + + CHK(s3d_scene_view_closest_point(scnview, pos, (float)INF, NULL, &hit) == RES_OK); + check_closest_point_cbox_sphere(pos, sphere_org, sphere_radius, walls_id, + sphere_id, NULL, 0, &hit); + } + + CHK(s3d_scene_view_ref_put(scnview) == RES_OK); + + /* Instantiate the cbox sphere scene */ + CHK(s3d_scene_instantiate(scn, &inst0) == RES_OK); + CHK(s3d_scene_instantiate(scn, &inst1) == RES_OK); + CHK(s3d_shape_get_id(inst0, &instances[0].id) == RES_OK); + CHK(s3d_shape_get_id(inst1, &instances[1].id) == RES_OK); + f3_mulf(instances[0].translation, sz, 0.5f); + CHK(s3d_instance_translate + (inst0, S3D_WORLD_TRANSFORM, instances[0].translation) == RES_OK); + f3_mulf(instances[1].translation, sz,-0.5f); + CHK(s3d_instance_translate + (inst1, S3D_WORLD_TRANSFORM, instances[1].translation) == RES_OK); + + /* Create a new scene with instantiated cbox sphere scenes */ + CHK(s3d_scene_ref_put(scn) == RES_OK); + CHK(s3d_scene_create(dev, &scn) == RES_OK); + CHK(s3d_scene_attach_shape(scn, inst0) == RES_OK); + CHK(s3d_scene_attach_shape(scn, inst1) == RES_OK); + + CHK(s3d_scene_view_create(scn, S3D_TRACE, &scnview) == RES_OK); + CHK(s3d_scene_view_get_aabb(scnview, low, upp) == RES_OK); + f3_mulf(mid, f3_add(mid, low, upp), 0.5f); + + /* Check point query on instances */ + FOR_EACH(i, 0, 10000) { + /* Randomly generate a point in a bounding box that is 2 times the size of + * the scene AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]); + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]); + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]); + + CHK(s3d_scene_view_closest_point(scnview, pos, (float)INF, NULL, &hit) == RES_OK); + check_closest_point_cbox_sphere(pos, sphere_org, sphere_radius, walls_id, + sphere_id, instances, 2/*#instances*/, &hit); + } + + /* Clean up */ + CHK(s3d_shape_ref_put(inst0) == RES_OK); + CHK(s3d_shape_ref_put(inst1) == RES_OK); + CHK(s3d_shape_ref_put(walls) == RES_OK); + CHK(s3d_shape_ref_put(sphere) == RES_OK); + CHK(s3d_scene_view_ref_put(scnview) == RES_OK); + CHK(s3d_scene_ref_put(scn) == RES_OK); +} + +/******************************************************************************* + * Cornell box test + ******************************************************************************/ +enum cbox_geom { + CBOX_WALLS, + CBOX_TALL_BLOCK, + CBOX_SHORT_BLOCK, + CBOX_GEOMS_COUNT__ +}; + +struct cbox_filter_data { + float query_pos[3]; + unsigned geom_to_filter[3]; +}; + +static int +cbox_filter + (const struct s3d_hit* hit, + const float org[3], + const float dir[3], + void* query_data, + void* filter_data) +{ + struct cbox_filter_data* data = query_data; + struct s3d_attrib attr; + float pos[3]; + float vec[3]; + + CHK(hit && org && dir && !S3D_HIT_NONE(hit)); + CHK((intptr_t)filter_data == (intptr_t)0xDECAFBAD); + CHK(f3_normalize(vec, dir) != 0); + + f3_add(pos, org, f3_mulf(pos, vec, hit->distance)); + CHK(s3d_primitive_get_attrib + (&hit->prim, S3D_POSITION, hit->uv, &attr) == RES_OK); + CHK(f3_eq_eps(attr.value, pos, POSITION_EPSILON)); + + if(!query_data) return 0; + + CHK(f3_eq_eps(data->query_pos, org, POSITION_EPSILON)); + + return data->geom_to_filter[0] == hit->prim.geom_id + || data->geom_to_filter[1] == hit->prim.geom_id + || data->geom_to_filter[2] == hit->prim.geom_id; +} + +static void +check_closest_point_cbox + (const float pos[3], + const unsigned geom_id[3], + struct s3d_hit* hit) +{ + struct closest_pt pt[CBOX_GEOMS_COUNT__] = { + CLOSEST_PT_NULL__, CLOSEST_PT_NULL__, CLOSEST_PT_NULL__ + }; + enum cbox_geom geom; + + CHK(pos && geom_id && hit); + + if(geom_id[CBOX_WALLS] != S3D_INVALID_ID) { /* Are the walls filtered */ + closest_point_mesh(pos, cbox_walls, cbox_walls_ids, cbox_walls_ntris, + geom_id[CBOX_WALLS], S3D_INVALID_ID, &pt[CBOX_WALLS]); + } + if(geom_id[CBOX_TALL_BLOCK] != S3D_INVALID_ID) { /* Is the block filtered */ + closest_point_mesh(pos, cbox_tall_block, cbox_block_ids, cbox_block_ntris, + geom_id[CBOX_TALL_BLOCK], S3D_INVALID_ID, &pt[CBOX_TALL_BLOCK]); + } + if(geom_id[CBOX_SHORT_BLOCK] != S3D_INVALID_ID) { /* Is the block filtered */ + closest_point_mesh(pos, cbox_short_block, cbox_block_ids, cbox_block_ntris, + geom_id[CBOX_SHORT_BLOCK], S3D_INVALID_ID, &pt[CBOX_SHORT_BLOCK]); + } + geom = pt[CBOX_WALLS].dst < pt[CBOX_TALL_BLOCK].dst + ? CBOX_WALLS : CBOX_TALL_BLOCK; + geom = pt[CBOX_SHORT_BLOCK].dst < pt[geom].dst + ? CBOX_SHORT_BLOCK : geom; + + if(pt[geom].dst >= FLT_MAX) { /* All geometries were filtered */ + CHK(S3D_HIT_NONE(hit)); + } else { + check_closest_point(hit, &pt[geom], 1); + } +} + +static void +test_cbox(struct s3d_device* dev) +{ + struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL; + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_scene* scn = NULL; + struct s3d_shape* walls = NULL; + struct s3d_shape* tall_block = NULL; + struct s3d_shape* short_block = NULL; + struct s3d_scene_view* scnview = NULL; + struct cbox_desc walls_desc; + struct cbox_desc tall_block_desc; + struct cbox_desc short_block_desc; + struct cbox_filter_data filter_data; + void* ptr = (void*)((intptr_t)0xDECAFBAD); + float pos[3]; + float low[3], upp[3], mid[3]; + unsigned geom_id[CBOX_GEOMS_COUNT__]; + size_t i; + + /* Create the Star-3D scene */ + CHK(s3d_scene_create(dev, &scn) == RES_OK); + CHK(s3d_shape_create_mesh(dev, &walls) == RES_OK); + CHK(s3d_shape_create_mesh(dev, &tall_block) == RES_OK); + CHK(s3d_shape_create_mesh(dev, &short_block) == RES_OK); + CHK(s3d_shape_get_id(walls, &geom_id[CBOX_WALLS]) == RES_OK); + CHK(s3d_shape_get_id(tall_block, &geom_id[CBOX_TALL_BLOCK]) == RES_OK); + CHK(s3d_shape_get_id(short_block, &geom_id[CBOX_SHORT_BLOCK]) == RES_OK); + CHK(s3d_mesh_set_hit_filter_function(walls, cbox_filter, ptr) == RES_OK); + CHK(s3d_mesh_set_hit_filter_function(tall_block, cbox_filter, ptr) == RES_OK); + CHK(s3d_mesh_set_hit_filter_function(short_block, cbox_filter, ptr) == RES_OK); + CHK(s3d_scene_attach_shape(scn, walls) == RES_OK); + CHK(s3d_scene_attach_shape(scn, tall_block) == RES_OK); + CHK(s3d_scene_attach_shape(scn, short_block) == RES_OK); + + vdata.usage = S3D_POSITION; + vdata.type = S3D_FLOAT3; + vdata.get = cbox_get_position; + + /* Setup the Cornell box walls */ + walls_desc.vertices = cbox_walls; + walls_desc.indices = cbox_walls_ids; + CHK(s3d_mesh_setup_indexed_vertices(walls, cbox_walls_ntris, cbox_get_ids, + cbox_walls_nverts, &vdata, 1, &walls_desc) == RES_OK); + + /* Setup the Cornell box tall block */ + tall_block_desc.vertices = cbox_tall_block; + tall_block_desc.indices = cbox_block_ids; + CHK(s3d_mesh_setup_indexed_vertices(tall_block, cbox_block_ntris, cbox_get_ids, + cbox_block_nverts, &vdata, 1, &tall_block_desc) == RES_OK); + + /* Setup the Cornell box short block */ + short_block_desc.vertices = cbox_short_block; + short_block_desc.indices = cbox_block_ids; + CHK(s3d_mesh_setup_indexed_vertices(short_block, cbox_block_ntris, cbox_get_ids, + cbox_block_nverts, &vdata, 1, &short_block_desc) == RES_OK); + + CHK(s3d_scene_view_create(scn, S3D_TRACE, &scnview) == RES_OK); + CHK(s3d_scene_view_get_aabb(scnview, low, upp) == RES_OK); + mid[0] = (low[0] + upp[0]) * 0.5f; + mid[1] = (low[1] + upp[1]) * 0.5f; + mid[2] = (low[2] + upp[2]) * 0.5f; + + /* Check point query on Cornell box */ + FOR_EACH(i, 0, 10000) { + /* Randomly generate a point in a bounding box that is 2 times the size of + * the triangle AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]); + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]); + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]); + + CHK(s3d_scene_view_closest_point(scnview, pos, (float)INF, NULL, &hit) == RES_OK); + check_closest_point_cbox(pos, geom_id, &hit); + } + + /* Filter the Cornell box blocks */ + filter_data.geom_to_filter[0] = geom_id[CBOX_TALL_BLOCK]; + filter_data.geom_to_filter[1] = geom_id[CBOX_SHORT_BLOCK]; + filter_data.geom_to_filter[2] = S3D_INVALID_ID; + geom_id[CBOX_TALL_BLOCK] = S3D_INVALID_ID; + geom_id[CBOX_SHORT_BLOCK] = S3D_INVALID_ID; + + /* Check point query filtering */ + FOR_EACH(i, 0, 10000) { + /* Randomly generate a point in a bounding box that is 2 times the size of + * the triangle AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]); + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]); + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]); + + f3_set(filter_data.query_pos, pos); + + CHK(s3d_scene_view_closest_point + (scnview, pos, (float)INF, &filter_data, &hit) == RES_OK); + + check_closest_point_cbox(pos, geom_id, &hit); + } + + /* Clean up */ + CHK(s3d_shape_ref_put(walls) == RES_OK); + CHK(s3d_shape_ref_put(tall_block) == RES_OK); + CHK(s3d_shape_ref_put(short_block) == RES_OK); + CHK(s3d_scene_ref_put(scn) == RES_OK); + CHK(s3d_scene_view_ref_put(scnview) == RES_OK); +} + +/******************************************************************************* + * Single triangle test + ******************************************************************************/ +static void +triangle_get_ids(const unsigned itri, unsigned ids[3], void* ctx) +{ + (void)ctx; + CHK(itri == 0); + CHK(ids); + ids[0] = 0; + ids[1] = 1; + ids[2] = 2; +} + +static void +triangle_get_pos(const unsigned ivert, float pos[3], void* ctx) +{ + (void)ctx; + CHK(ivert < 3); + CHK(pos); + switch(ivert) { /* Setup a random triangle */ + case 0: f3(pos, -0.5f, -0.3f, 0.1f); break; + case 1: f3(pos, -0.4f, 0.2f, 0.3f); break; + case 2: f3(pos, 0.7f, 0.01f, -0.5f); break; + default: FATAL("Unreachable code\n"); break; + } +} + +static void +test_single_triangle(struct s3d_device* dev) +{ + struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL; + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_scene* scn = NULL; + struct s3d_scene_view* view = NULL; + struct s3d_shape* msh = NULL; + struct s3d_attrib attr; + float v0[3], v1[3], v2[3]; + float pos[3] = {0,0,0}; + float closest_pos[3] = {0,0,0}; + float low[3], upp[3], mid[3]; + size_t i; + + CHK(s3d_scene_create(dev, &scn) == RES_OK); + CHK(s3d_shape_create_mesh(dev, &msh) == RES_OK); + CHK(s3d_scene_attach_shape(scn, msh) == RES_OK); + + vdata.usage = S3D_POSITION; + vdata.type = S3D_FLOAT3; + vdata.get = triangle_get_pos; + CHK(s3d_mesh_setup_indexed_vertices + (msh, 1, triangle_get_ids, 3, &vdata, 1, NULL) == RES_OK); + + CHK(s3d_scene_view_create(scn, S3D_TRACE, &view) == RES_OK); + + triangle_get_pos(0, v0, NULL); + triangle_get_pos(1, v1, NULL); + triangle_get_pos(2, v2, NULL); + + /* Compute the triangle AABB */ + low[0] = MMIN(MMIN(v0[0], v1[0]), v2[0]); + low[1] = MMIN(MMIN(v0[1], v1[1]), v2[1]); + low[2] = MMIN(MMIN(v0[2], v1[2]), v2[2]); + upp[0] = MMAX(MMAX(v0[0], v1[0]), v2[0]); + upp[1] = MMAX(MMAX(v0[1], v1[1]), v2[1]); + upp[2] = MMAX(MMAX(v0[2], v1[2]), v2[2]); + mid[0] = (low[0] + upp[0]) * 0.5f; + mid[1] = (low[1] + upp[1]) * 0.5f; + mid[2] = (low[2] + upp[2]) * 0.5f; + + FOR_EACH(i, 0, 10000) { + /* Randomly generate a point in a bounding box that is 10 times the size of + * the triangle AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]) * 5.f; + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]) * 5.f; + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]) * 5.f; + + CHK(s3d_scene_view_closest_point(view, pos, (float)INF, NULL, &hit) == RES_OK); + CHK(!S3D_HIT_NONE(&hit)); + CHK(s3d_primitive_get_attrib(&hit.prim, S3D_POSITION, hit.uv, &attr) == RES_OK); + + /* Cross check the point query result */ + closest_point_triangle(pos, v0, v1, v2, closest_pos); + CHK(f3_eq_eps(closest_pos, attr.value, 1.e-4f)); + } + + FOR_EACH(i, 0, 10000) { + float radius; + + /* Randomly generate a point in a bounding box that is 10 times the size of + * the triangle AABB */ + pos[0] = mid[0] + (rand_canonic() * 2 - 1) * (upp[0] - low[0]) * 5.f; + pos[1] = mid[1] + (rand_canonic() * 2 - 1) * (upp[1] - low[1]) * 5.f; + pos[2] = mid[2] + (rand_canonic() * 2 - 1) * (upp[2] - low[2]) * 5.f; + + CHK(s3d_scene_view_closest_point(view, pos, (float)INF, NULL, &hit) == RES_OK); + CHK(!S3D_HIT_NONE(&hit)); + + /* Check that the radius is an exclusive upper bound */ + radius = hit.distance; + CHK(s3d_scene_view_closest_point(view, pos, radius, NULL, &hit) == RES_OK); + CHK(S3D_HIT_NONE(&hit)); + radius = nextafterf(radius, FLT_MAX); + CHK(s3d_scene_view_closest_point(view, pos, radius, NULL, &hit) == RES_OK); + CHK(!S3D_HIT_NONE(&hit)); + CHK(hit.distance == nextafterf(radius, 0.f)); + } + + CHK(s3d_shape_ref_put(msh) == RES_OK); + CHK(s3d_scene_view_ref_put(view) == RES_OK); + CHK(s3d_scene_ref_put(scn) == RES_OK); +} + +/******************************************************************************* + * Miscellaneous test + ******************************************************************************/ +static void +test_api(struct s3d_device* dev) +{ + struct s3d_hit hit = S3D_HIT_NULL; + struct s3d_scene* scn = NULL; + struct s3d_scene_view* view = NULL; + float pos[3] = {0,0,0}; + + CHK(s3d_scene_create(dev, &scn) == RES_OK); + CHK(s3d_scene_view_create(scn, S3D_TRACE, &view) == RES_OK); + + CHK(s3d_scene_view_closest_point(NULL, pos, 1.f, NULL, &hit) == RES_BAD_ARG); + CHK(s3d_scene_view_closest_point(view, NULL, 1.f, NULL, &hit) == RES_BAD_ARG); + CHK(s3d_scene_view_closest_point(view, pos, 0.f, NULL, &hit) == RES_BAD_ARG); + CHK(s3d_scene_view_closest_point(view, pos, 1.f, NULL, NULL) == RES_BAD_ARG); + CHK(s3d_scene_view_closest_point(view, pos, 1.f, NULL, &hit) == RES_OK); + CHK(S3D_HIT_NONE(&hit)); + + CHK(s3d_scene_view_ref_put(view) == RES_OK); + CHK(s3d_scene_view_create(scn, S3D_SAMPLE, &view) == RES_OK); + CHK(s3d_scene_view_closest_point(view, pos, 1.f, NULL, &hit) == RES_BAD_OP); + + CHK(s3d_scene_view_ref_put(view) == RES_OK); + CHK(s3d_scene_ref_put(scn) == RES_OK); +} + +/******************************************************************************* + * Main function + ******************************************************************************/ +int +main(int argc, char** argv) +{ + struct mem_allocator allocator; + struct s3d_device* dev = NULL; + (void)argc, (void)argv; + + mem_init_proxy_allocator(&allocator, &mem_default_allocator); + CHK(s3d_device_create(NULL, &allocator, 1, &dev) == RES_OK); + + test_api(dev); + test_single_triangle(dev); + test_cbox(dev); + test_cbox_sphere(dev); + + CHK(s3d_device_ref_put(dev) == RES_OK); + + check_memory_allocator(&allocator); + mem_shutdown_proxy_allocator(&allocator); + CHK(mem_allocated_size() == 0); + return 0; +}