htrdr

Solving radiative transfer in heterogeneous media
git clone git://git.meso-star.fr/htrdr.git
Log | Files | Refs | README | LICENSE

htrdr_geometry.c (21279B)


      1 /* Copyright (C) 2018-2019, 2022-2025 Centre National de la Recherche Scientifique
      2  * Copyright (C) 2020-2022 Institut Mines Télécom Albi-Carmaux
      3  * Copyright (C) 2022-2025 Institut Pierre-Simon Laplace
      4  * Copyright (C) 2022-2025 Institut de Physique du Globe de Paris
      5  * Copyright (C) 2018-2025 |Méso|Star> (contact@meso-star.com)
      6  * Copyright (C) 2022-2025 Observatoire de Paris
      7  * Copyright (C) 2022-2025 Université de Reims Champagne-Ardenne
      8  * Copyright (C) 2022-2025 Université de Versaille Saint-Quentin
      9  * Copyright (C) 2018-2019, 2022-2025 Université Paul Sabatier
     10  *
     11  * This program is free software: you can redistribute it and/or modify
     12  * it under the terms of the GNU General Public License as published by
     13  * the Free Software Foundation, either version 3 of the License, or
     14  * (at your option) any later version.
     15  *
     16  * This program is distributed in the hope that it will be useful,
     17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     19  * GNU General Public License for more details.
     20  *
     21  * You should have received a copy of the GNU General Public License
     22  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     23 
     24 #define _POSIX_C_SOURCE 200112L /* strtok_r support */
     25 
     26 #include "core/htrdr.h"
     27 #include "core/htrdr_geometry.h"
     28 #include "core/htrdr_interface.h"
     29 #include "core/htrdr_log.h"
     30 #include "core/htrdr_materials.h"
     31 #include "core/htrdr_slab.h"
     32 
     33 #include <aw.h>
     34 #include <star/s3d.h>
     35 
     36 #include <rsys/clock_time.h>
     37 #include <rsys/cstr.h>
     38 #include <rsys/dynamic_array_double.h>
     39 #include <rsys/dynamic_array_size_t.h>
     40 #include <rsys/double2.h>
     41 #include <rsys/double3.h>
     42 #include <rsys/float2.h>
     43 #include <rsys/float3.h>
     44 #include <rsys/hash_table.h>
     45 #include <rsys/ref_count.h>
     46 
     47 #include <string.h> /* strtok_r */
     48 
     49 /* Define the hash table that maps an Obj vertex id to its position into the
     50  * vertex buffer */
     51 #define HTABLE_NAME vertex
     52 #define HTABLE_KEY size_t /* Obj vertex id */
     53 #define HTABLE_DATA size_t
     54 #include <rsys/hash_table.h>
     55 
     56 /* Define the hash table that maps the Star-3D shape id to its interface */
     57 #define HTABLE_NAME interface
     58 #define HTABLE_KEY unsigned /* Star-3D shape id */
     59 #define HTABLE_DATA struct htrdr_interface
     60 #include <rsys/hash_table.h>
     61 
     62 struct mesh {
     63   const struct darray_double* positions;
     64   const struct darray_size_t* indices;
     65 };
     66 static const struct mesh MESH_NULL;
     67 
     68 struct ray_context {
     69   struct s3d_hit hit_from;
     70 
     71   s3d_hit_filter_function_T user_filter; /* May be NULL */
     72   void* user_filter_data; /* May be NULL */
     73 };
     74 #define RAY_CONTEXT_NULL__ {S3D_HIT_NULL__, NULL, NULL}
     75 static const struct ray_context RAY_CONTEXT_NULL = RAY_CONTEXT_NULL__;
     76 
     77 struct htrdr_geometry {
     78   struct s3d_scene_view* view;
     79   float lower[3]; /* Ground lower bound */
     80   float upper[3]; /* Ground upper bound */
     81   int repeat; /* Make the geom infinite in X and Y */
     82 
     83   /* A empirical value relative to the extent of the geometry that represents
     84    * the threshold below which a numerical problem could occur. */
     85   float epsilon;
     86 
     87   struct htable_interface interfaces; /* Map a Star3D shape to its interface */
     88 
     89   struct htrdr* htrdr;
     90   ref_T ref;
     91 };
     92 
     93 /*******************************************************************************
     94  * Helper functions
     95  ******************************************************************************/
     96 /* Check that `hit' roughly lies on an edge. For triangular primitives, a
     97  * simple but approximative way is to test that its position have at least one
     98  * barycentric coordinate roughly equal to 0 or 1. */
     99 static FINLINE int
    100 hit_on_edge(const struct s3d_hit* hit)
    101 {
    102   const float on_edge_eps = 1.e-4f;
    103   float w;
    104   ASSERT(hit && !S3D_HIT_NONE(hit));
    105   w = 1.f - hit->uv[0] - hit->uv[1];
    106   return eq_epsf(hit->uv[0], 0.f, on_edge_eps)
    107       || eq_epsf(hit->uv[0], 1.f, on_edge_eps)
    108       || eq_epsf(hit->uv[1], 0.f, on_edge_eps)
    109       || eq_epsf(hit->uv[1], 1.f, on_edge_eps)
    110       || eq_epsf(w, 0.f, on_edge_eps)
    111       || eq_epsf(w, 1.f, on_edge_eps);
    112 }
    113 
    114 /* Return 1 if the submitted hit is actually a self hit, i.e. if the ray starts
    115  * the primitive from which it starts */
    116 static INLINE int
    117 self_hit
    118   (const struct s3d_hit* hit,
    119    const float ray_org[3],
    120    const float ray_dir[3],
    121    const float ray_range[2],
    122    const struct s3d_hit* hit_from)
    123 {
    124   ASSERT(hit && hit_from);
    125   (void)ray_org, (void)ray_dir;
    126 
    127   /* Internally, Star-3D relies on Embree that, due to numerically inaccuracy,
    128    * can find hits whose distances are not strictly included in the submitted
    129    * ray range. Discard these hits. */
    130   if(hit->distance <= ray_range[0] || hit->distance >= ray_range[1])
    131     return 1;
    132 
    133   /* The hit primitive is the one from which the the ray starts. Discard these
    134    * hits */
    135   if(S3D_PRIMITIVE_EQ(&hit->prim, &hit_from->prim))
    136     return 1;
    137 
    138   /* If the targeted point is near of the origin, check that it lies on an
    139    * edge/vertex shared by the 2 primitives. In such situation, we assume that
    140    * it is a self intersection and discard this hit */
    141   if(!S3D_HIT_NONE(hit_from)
    142   && eq_epsf(hit->distance, 0, 1.e-1f)
    143   && hit_on_edge(hit_from) 
    144   && hit_on_edge(hit)) {
    145     return 1;
    146   }
    147 
    148   /* No self hit */
    149   return 0;
    150 }
    151 
    152 static int
    153 geometry_filter
    154   (const struct s3d_hit* hit,
    155    const float ray_org[3],
    156    const float ray_dir[3],
    157    const float ray_range[2],
    158    void* ray_data,
    159    void* filter_data)
    160 {
    161   const struct ray_context* ray_ctx = ray_data;
    162   (void)filter_data;
    163 
    164   if(!ray_ctx) /* Nothing to do */
    165     return 0;
    166 
    167   if(self_hit(hit, ray_org, ray_dir, ray_range, &ray_ctx->hit_from))
    168     return 1; /* Discard this hit */
    169 
    170   if(!ray_ctx->user_filter) /* That's all */
    171     return 0;
    172 
    173   return ray_ctx->user_filter /* Invoke user filtering */
    174     (hit, ray_org, ray_dir, ray_range, ray_ctx->user_filter_data, filter_data);
    175 }
    176 
    177 static res_T
    178 parse_shape_interface
    179   (struct htrdr* htrdr,
    180    struct htrdr_materials* mats,
    181    const char* name,
    182    struct htrdr_interface* interf)
    183 {
    184   struct str str;
    185   char* mtl_name0 = NULL;
    186   char* mtl_name1 = NULL;
    187   char* mtl_name2 = NULL;
    188   char* mtl_name_front = NULL;
    189   char* mtl_name_thin = NULL;
    190   char* mtl_name_back = NULL;
    191   char* tk_ctx = NULL;
    192   int has_front = 0;
    193   int has_thin = 0;
    194   int has_back = 0;
    195   res_T res = RES_OK;
    196   ASSERT(htrdr && mats && name && interf);
    197 
    198   str_init(htrdr_get_allocator(htrdr), &str);
    199 
    200   /* Locally copy the string to parse */
    201   res = str_set(&str, name);
    202   if(res != RES_OK) {
    203     htrdr_log_err(htrdr,
    204       "Could not locally copy the shape material string `%s' -- %s.\n",
    205       name, res_to_cstr(res));
    206     goto error;
    207   }
    208 
    209   /* Reset the interface */
    210   memset(interf, 0, sizeof(*interf));
    211 
    212   /* Parse the name of the front/back/thin materials */
    213   mtl_name0 = strtok_r(str_get(&str), ":", &tk_ctx);
    214   mtl_name1 = strtok_r(NULL, ":", &tk_ctx);
    215   mtl_name2 = strtok_r(NULL, ":", &tk_ctx);
    216   ASSERT(mtl_name0); /* This can't be NULL */
    217   mtl_name_front = mtl_name0;
    218   if(mtl_name2) {
    219     mtl_name_thin = mtl_name1;
    220     mtl_name_back = mtl_name2;
    221   } else {
    222     mtl_name_thin = NULL;
    223     mtl_name_back = mtl_name1;
    224   }
    225 
    226   if(!mtl_name_back) {
    227     htrdr_log_err(htrdr,
    228       "The back material name is missing `%s'.\n", name);
    229     res = RES_BAD_ARG;
    230     goto error;
    231   }
    232 
    233   /* Fetch the interface material if any */
    234   if(mtl_name_thin) {
    235     has_thin = htrdr_materials_find_mtl(mats, mtl_name_thin, &interf->mtl_thin);
    236     if(!has_thin) {
    237       htrdr_log_err(htrdr,
    238         "Invalid interface `%s'. The interface material `%s' is unknown.\n",
    239         name, mtl_name_thin);
    240       res = RES_BAD_ARG;
    241       goto error;
    242     }
    243   }
    244 
    245   /* Fetch the front material */
    246   has_front = htrdr_materials_find_mtl(mats, mtl_name_front, &interf->mtl_front);
    247   if(!has_front) {
    248     htrdr_log_err(htrdr,
    249       "Invalid interface `%s'. The front material `%s' is unknown.\n",
    250       name, mtl_name_front);
    251     res = RES_BAD_ARG;
    252     goto error;
    253   }
    254 
    255   /* Fetch the back material */
    256   has_back = htrdr_materials_find_mtl(mats, mtl_name_back, &interf->mtl_back);
    257   if(!has_back) {
    258     htrdr_log_err(htrdr,
    259       "Invalid interface `%s'. The back material `%s' is unknown.\n",
    260       name, mtl_name_back);
    261     res = RES_BAD_ARG;
    262     goto error;
    263   }
    264 
    265 exit:
    266   str_release(&str);
    267   return res;
    268 error:
    269   *interf = HTRDR_INTERFACE_NULL;
    270   goto exit;
    271 }
    272 
    273 static res_T
    274 setup_mesh
    275   (struct htrdr* htrdr,
    276    const char* filename,
    277    struct aw_obj* obj,
    278    struct aw_obj_named_group* mtl,
    279    struct darray_double* positions,
    280    struct darray_size_t* indices,
    281    struct htable_vertex* vertices) /* Scratch data structure */
    282 {
    283   size_t iface;
    284   res_T res = RES_OK;
    285   ASSERT(htrdr && filename && obj && mtl && positions && indices && vertices);
    286 
    287   darray_double_clear(positions);
    288   darray_size_t_clear(indices);
    289   htable_vertex_clear(vertices);
    290 
    291   FOR_EACH(iface, mtl->face_id, mtl->face_id+mtl->faces_count) {
    292     struct aw_obj_face face;
    293     size_t ivertex;
    294 
    295     AW(obj_get_face(obj, iface, &face));
    296     if(face.vertices_count != 3) {
    297       htrdr_log_err(htrdr,
    298         "The obj `%s' has non-triangulated polygons "
    299         "while currently only triangular meshes are supported.\n",
    300         filename);
    301       res = RES_BAD_ARG;
    302       goto error;
    303     }
    304 
    305     FOR_EACH(ivertex, 0, face.vertices_count) {
    306       struct aw_obj_vertex vertex;
    307       size_t id;
    308       size_t* pid;
    309 
    310       AW(obj_get_vertex(obj, face.vertex_id + ivertex, &vertex));
    311       pid = htable_vertex_find(vertices, &vertex.position_id);
    312       if(pid) {
    313         id = *pid;
    314       } else {
    315         struct aw_obj_vertex_data vdata;
    316 
    317         id = darray_double_size_get(positions) / 3;
    318 
    319         res = darray_double_resize(positions, id*3 + 3);
    320         if(res != RES_OK) {
    321           htrdr_log_err(htrdr,
    322             "Could not locally copy the vertex position -- %s.\n",
    323             res_to_cstr(res));
    324           goto error;
    325         }
    326 
    327         AW(obj_get_vertex_data(obj, &vertex, &vdata));
    328         darray_double_data_get(positions)[id*3+0] = vdata.position[0];
    329         darray_double_data_get(positions)[id*3+1] = vdata.position[1];
    330         darray_double_data_get(positions)[id*3+2] = vdata.position[2];
    331 
    332         res = htable_vertex_set(vertices, &vertex.position_id, &id);
    333         if(res != RES_OK) {
    334           htrdr_log_err(htrdr,
    335             "Could not register the vertex position -- %s.\n",
    336              res_to_cstr(res));
    337           goto error;
    338         }
    339       }
    340 
    341       res = darray_size_t_push_back(indices, &id);
    342       if(res != RES_OK) {
    343         htrdr_log_err(htrdr,
    344           "Could not locally copy the face index -- %s\n",
    345           res_to_cstr(res));
    346         goto error;
    347       }
    348     }
    349   }
    350 exit:
    351   return res;
    352 error:
    353   darray_double_clear(positions);
    354   darray_size_t_clear(indices);
    355   htable_vertex_clear(vertices);
    356   goto exit;
    357 }
    358 
    359 static void
    360 get_position(const unsigned ivert, float position[3], void* ctx)
    361 {
    362   const struct mesh* mesh = ctx;
    363   const double* pos = NULL;
    364   ASSERT(mesh);
    365   ASSERT(ivert < darray_double_size_get(mesh->positions) / 3);
    366   pos = darray_double_cdata_get(mesh->positions) + ivert*3;
    367   position[0] = (float)pos[0];
    368   position[1] = (float)pos[1];
    369   position[2] = (float)pos[2];
    370 }
    371 
    372 static void
    373 get_indices(const unsigned itri, unsigned indices[3], void* ctx)
    374 {
    375   const struct mesh* mesh = ctx;
    376   const size_t* ids = NULL;
    377   ASSERT(mesh);
    378   ASSERT(itri < darray_size_t_size_get(mesh->indices) / 3);
    379   ids = darray_size_t_cdata_get(mesh->indices) + itri*3;
    380   indices[0] = (unsigned)ids[0];
    381   indices[1] = (unsigned)ids[1];
    382   indices[2] = (unsigned)ids[2];
    383 }
    384 
    385 static res_T
    386 create_s3d_shape
    387   (struct htrdr* htrdr,
    388    const struct darray_double* positions,
    389    const struct darray_size_t* indices,
    390    struct s3d_shape** out_shape)
    391 {
    392   struct s3d_shape* shape = NULL;
    393   struct s3d_vertex_data vdata = S3D_VERTEX_DATA_NULL;
    394   struct mesh mesh = MESH_NULL;
    395   res_T res = RES_OK;
    396   ASSERT(htrdr && positions && indices && out_shape);
    397   ASSERT(darray_double_size_get(positions) != 0);
    398   ASSERT(darray_size_t_size_get(indices) != 0);
    399   ASSERT(darray_double_size_get(positions)%3 == 0);
    400   ASSERT(darray_size_t_size_get(indices)%3 == 0);
    401 
    402   mesh.positions = positions;
    403   mesh.indices = indices;
    404 
    405   res = s3d_shape_create_mesh(htrdr_get_s3d(htrdr), &shape);
    406   if(res != RES_OK) {
    407     htrdr_log_err(htrdr, "Error creating a Star-3D shape -- %s.\n",
    408       res_to_cstr(res));
    409     goto error;
    410   }
    411 
    412   vdata.usage = S3D_POSITION;
    413   vdata.type = S3D_FLOAT3;
    414   vdata.get = get_position;
    415 
    416   res = s3d_mesh_setup_indexed_vertices
    417     (shape,
    418      (unsigned int)(darray_size_t_size_get(indices)/3),
    419      get_indices,
    420      (unsigned int)(darray_double_size_get(positions)/3),
    421      &vdata, 1,
    422      &mesh);
    423   if(res != RES_OK){
    424     htrdr_log_err(htrdr, "Could not setup the Star-3D shape -- %s.\n",
    425       res_to_cstr(res));
    426     goto error;
    427   }
    428 
    429   res = s3d_mesh_set_hit_filter_function(shape, geometry_filter, NULL);
    430   if(res != RES_OK) {
    431     htrdr_log_err(htrdr,
    432       "Could not setup the Star-3D hit filter function of the geometry "
    433        "-- %s.\n", res_to_cstr(res));
    434     goto error;
    435   }
    436 
    437 exit:
    438   *out_shape = shape;
    439   return res;
    440 error:
    441   if(shape) {
    442     S3D(shape_ref_put(shape));
    443     shape = NULL;
    444   }
    445   goto exit;
    446 }
    447 
    448 static res_T
    449 setup_geometry
    450   (struct htrdr_geometry* geom,
    451    struct htrdr_materials* mats,
    452    const char* obj_filename)
    453 {
    454   struct aw_obj_desc desc;
    455   struct htable_vertex vertices;
    456   struct darray_double positions;
    457   struct darray_size_t indices;
    458   struct aw_obj* obj = NULL;
    459   struct s3d_shape* shape = NULL;
    460   struct s3d_scene* scene = NULL;
    461   size_t iusemtl;
    462 
    463   res_T res = RES_OK;
    464   ASSERT(geom && mats && obj_filename);
    465 
    466   htable_vertex_init(htrdr_get_allocator(geom->htrdr), &vertices);
    467   darray_double_init(htrdr_get_allocator(geom->htrdr), &positions);
    468   darray_size_t_init(htrdr_get_allocator(geom->htrdr), &indices);
    469 
    470   res = aw_obj_create
    471     (htrdr_get_logger(geom->htrdr),
    472      htrdr_get_allocator(geom->htrdr),
    473      htrdr_get_verbosity_level(geom->htrdr),
    474      &obj);
    475   if(res != RES_OK) {
    476     htrdr_log_err(geom->htrdr, "Could not create the obj loader -- %s.\n",
    477       res_to_cstr(res));
    478     goto error;
    479   }
    480 
    481   res = s3d_scene_create(htrdr_get_s3d(geom->htrdr), &scene);
    482   if(res != RES_OK) {
    483     htrdr_log_err(geom->htrdr, "Could not create the Star-3D scene -- %s.\n",
    484       res_to_cstr(res));
    485     goto error;
    486   }
    487 
    488   /* Load the geometry data */
    489   res = aw_obj_load(obj, obj_filename);
    490   if(res != RES_OK) goto error;
    491 
    492   /* Fetch the descriptor of the loaded geometry */
    493   AW(obj_get_desc(obj, &desc));
    494 
    495   if(desc.usemtls_count == 0) {
    496     htrdr_log_err(geom->htrdr, "The obj `%s' has no material.\n", obj_filename);
    497     res = RES_BAD_ARG;
    498     goto error;
    499   }
    500 
    501   /* Setup the geometry */
    502   FOR_EACH(iusemtl, 0, desc.usemtls_count) {
    503     struct aw_obj_named_group mtl;
    504     struct htrdr_interface interf;
    505     unsigned ishape;
    506 
    507     AW(obj_get_mtl(obj, iusemtl , &mtl));
    508 
    509     res = parse_shape_interface(geom->htrdr, mats, mtl.name, &interf);
    510     if(res != RES_OK) goto error;
    511 
    512     res = setup_mesh
    513       (geom->htrdr, obj_filename, obj, &mtl, &positions, &indices, &vertices);
    514     if(res != RES_OK) goto error;
    515 
    516     res = create_s3d_shape(geom->htrdr, &positions, &indices, &shape);
    517     if(res != RES_OK) goto error;
    518 
    519     S3D(shape_get_id(shape, &ishape));
    520     res = htable_interface_set(&geom->interfaces, &ishape, &interf);
    521     if(res != RES_OK) {
    522       htrdr_log_err(geom->htrdr,
    523         "Could not map the Star-3D shape to its Star-Materials -- %s.\n",
    524         res_to_cstr(res));
    525       goto error;
    526     }
    527 
    528     res = s3d_scene_attach_shape(scene, shape);
    529     if(res != RES_OK) {
    530       htrdr_log_err(geom->htrdr,
    531         "Could not attach a Star-3D shape to the Star-3D scene -- %s.\n",
    532         res_to_cstr(res));
    533       goto error;
    534     }
    535 
    536     S3D(shape_ref_put(shape));
    537     shape = NULL;
    538   }
    539 
    540   res = s3d_scene_view_create(scene, S3D_GET_PRIMITIVE|S3D_TRACE, &geom->view);
    541   if(res != RES_OK) {
    542     htrdr_log_err(geom->htrdr,
    543       "Could not create the Star-3D scene view -- %s.\n",
    544       res_to_cstr(res));
    545     goto error;
    546   }
    547 
    548   res = s3d_scene_view_get_aabb(geom->view, geom->lower, geom->upper);
    549   if(res != RES_OK) {
    550     htrdr_log_err(geom->htrdr,
    551       "Could not get the bounding box of the geometry -- %s.\n",
    552       res_to_cstr(res));
    553     goto error;
    554   }
    555 
    556 exit:
    557   if(obj) AW(obj_ref_put(obj));
    558   if(shape) S3D(shape_ref_put(shape));
    559   if(scene) S3D(scene_ref_put(scene));
    560   htable_vertex_release(&vertices);
    561   darray_double_release(&positions);
    562   darray_size_t_release(&indices);
    563   return res;
    564 error:
    565   goto exit;
    566 }
    567 
    568 static void
    569 release_ground(ref_T* ref)
    570 {
    571   struct htrdr_geometry* geom;
    572   struct htrdr* htrdr;
    573   ASSERT(ref);
    574   geom = CONTAINER_OF(ref, struct htrdr_geometry, ref);
    575   if(geom->view) S3D(scene_view_ref_put(geom->view));
    576   htable_interface_release(&geom->interfaces);
    577   htrdr = geom->htrdr;
    578   MEM_RM(htrdr_get_allocator(geom->htrdr), geom);
    579   htrdr_ref_put(htrdr);
    580 }
    581 
    582 /*******************************************************************************
    583  * Local functions
    584  ******************************************************************************/
    585 res_T
    586 htrdr_geometry_create
    587   (struct htrdr* htrdr,
    588    const char* obj_filename, /* May be NULL */
    589    struct htrdr_materials* mats,
    590    struct htrdr_geometry** out_ground)
    591 {
    592   char buf[128];
    593   struct htrdr_geometry* geom = NULL;
    594   double low[3];
    595   double upp[3];
    596   double tmp[3];
    597   double extent;
    598   struct time t0, t1;
    599   res_T res = RES_OK;
    600   ASSERT(htrdr && obj_filename && mats && out_ground);
    601 
    602   geom = MEM_CALLOC(htrdr_get_allocator(htrdr), 1, sizeof(*geom));
    603   if(!geom) {
    604     res = RES_MEM_ERR;
    605     htrdr_log_err(htrdr,
    606       "%s: could not allocate the geom data structure -- %s.\n",
    607       FUNC_NAME, res_to_cstr(res));
    608     goto error;
    609   }
    610   ref_init(&geom->ref);
    611   f3_splat(geom->lower, (float)INF);
    612   f3_splat(geom->upper,-(float)INF);
    613   htable_interface_init(htrdr_get_allocator(htrdr), &geom->interfaces);
    614   htrdr_ref_get(htrdr);
    615   geom->htrdr = htrdr;
    616 
    617   htrdr_log(geom->htrdr, "Loading geometry from `%s'.\n", obj_filename);
    618   time_current(&t0);
    619   res = setup_geometry(geom, mats, obj_filename);
    620   if(res != RES_OK) goto error;
    621   time_sub(&t0, time_current(&t1), &t0);
    622   time_dump(&t0, TIME_ALL, NULL, buf, sizeof(buf));
    623   htrdr_log(geom->htrdr, "Setup geom in %s\n", buf);
    624 
    625   htrdr_geometry_get_aabb(geom, low, upp);
    626   extent = d3_len(d3_sub(tmp, upp, low));
    627   geom->epsilon = MMAX((float)(extent * 1e-6), FLT_EPSILON);
    628 
    629 exit:
    630   *out_ground = geom;
    631   return res;
    632 error:
    633   if(geom) {
    634     htrdr_geometry_ref_put(geom);
    635     geom = NULL;
    636   }
    637   goto exit;
    638 }
    639 
    640 void
    641 htrdr_geometry_ref_get(struct htrdr_geometry* geom)
    642 {
    643   ASSERT(geom);
    644   ref_get(&geom->ref);
    645 }
    646 
    647 void
    648 htrdr_geometry_ref_put(struct htrdr_geometry* geom)
    649 {
    650   ASSERT(geom);
    651   ref_put(&geom->ref, release_ground);
    652 }
    653 
    654 void
    655 htrdr_geometry_get_interface
    656   (struct htrdr_geometry* geom,
    657    const struct s3d_hit* hit,
    658    struct htrdr_interface* out_interface)
    659 {
    660   struct htrdr_interface* interf = NULL;
    661   ASSERT(geom && hit && out_interface);
    662 
    663   interf = htable_interface_find(&geom->interfaces, &hit->prim.geom_id);
    664   ASSERT(interf);
    665 
    666   *out_interface = *interf;
    667 }
    668 
    669 void
    670 htrdr_geometry_get_hit_position
    671   (const struct htrdr_geometry* geom,
    672    const struct s3d_hit* hit,
    673    double position[3])
    674 {
    675   struct s3d_attrib attr;
    676   ASSERT(geom && hit && position && !S3D_HIT_NONE(hit));
    677   (void)geom;
    678 
    679   S3D(primitive_get_attrib(&hit->prim, S3D_POSITION, hit->uv, &attr));
    680   position[0] = attr.value[0];
    681   position[1] = attr.value[1];
    682   position[2] = attr.value[2];
    683 }
    684 
    685 res_T
    686 htrdr_geometry_trace_ray
    687   (struct htrdr_geometry* geom,
    688    const struct htrdr_geometry_trace_ray_args* args,
    689    struct s3d_hit* hit)
    690 {
    691   struct ray_context ray_ctx = RAY_CONTEXT_NULL;
    692   float ray_org[3];
    693   float ray_dir[3];
    694   float ray_range[2];
    695   res_T res = RES_OK;
    696   ASSERT(geom && args && hit);
    697 
    698   f3_set_d3(ray_org, args->ray_org);
    699   f3_set_d3(ray_dir, args->ray_dir);
    700   f2_set_d2(ray_range, args->ray_range);
    701   ray_ctx.hit_from = args->hit_from;
    702   ray_ctx.user_filter = args->filter;
    703   ray_ctx.user_filter_data = args->filter_context;
    704 
    705   res = s3d_scene_view_trace_ray
    706     (geom->view, ray_org, ray_dir, ray_range, &ray_ctx, hit);
    707   if(res != RES_OK) {
    708     htrdr_log_err(geom->htrdr,
    709       "%s: could not trace the ray against the geometry -- %s.\n",
    710       FUNC_NAME, res_to_cstr(res));
    711     goto error;
    712   }
    713 
    714 exit:
    715   return res;
    716 error:
    717   goto exit;
    718 }
    719 
    720 res_T
    721 htrdr_geometry_find_closest_point
    722   (struct htrdr_geometry* geom,
    723    const double pos[3],
    724    const double radius,
    725    struct s3d_hit* hit)
    726 {
    727   float query_pos[3];
    728   float query_radius;
    729   res_T res = RES_OK;
    730   ASSERT(geom && pos && hit);
    731 
    732   query_radius = (float)radius;
    733   f3_set_d3(query_pos, pos);
    734 
    735   /* Closest point query */
    736   res = s3d_scene_view_closest_point
    737     (geom->view, query_pos, query_radius, NULL, hit);
    738   if(res != RES_OK) {
    739     htrdr_log_err(geom->htrdr,
    740       "%s: could not query the closest point to the geometry -- %s.\n",
    741       FUNC_NAME, res_to_cstr(res));
    742     goto error;
    743   }
    744 
    745 exit:
    746   return res;
    747 error:
    748   goto exit;
    749 }
    750 
    751 void
    752 htrdr_geometry_get_aabb
    753   (const struct htrdr_geometry* geom,
    754    double lower[3],
    755    double upper[3])
    756 {
    757   ASSERT(geom && lower && upper);
    758   d3_set_f3(lower, geom->lower);
    759   d3_set_f3(upper, geom->upper);
    760 }
    761 
    762 double
    763 htrdr_geometry_get_epsilon(const struct htrdr_geometry* geom)
    764 {
    765   ASSERT(geom);
    766   return geom->epsilon;
    767 }