star-cad

Geometric operators for computer-aided design
git clone git://git.meso-star.fr/star-cad.git
Log | Files | Refs | README | LICENSE

scad.c (26955B)


      1 /* Copyright (C) 2022-2024 |Méso|Star> (contact@meso-star.com)
      2  *
      3  * This program is free software: you can redistribute it and/or modify
      4  * it under the terms of the GNU General Public License as published by
      5  * the Free Software Foundation, either version 3 of the License, or
      6  * (at your option) any later version.
      7  *
      8  * This program is distributed in the hope that it will be useful,
      9  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     11  * GNU General Public License for more details.
     12  *
     13  * You should have received a copy of the GNU General Public License
     14  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     15 
     16 #include "scad.h"
     17 #include "scad_c.h"
     18 #include "scad_device.h"
     19 #include "scad_geometry.h"
     20 
     21 #include <gmsh/gmshc.h>
     22 
     23 #include <rsys/cstr.h>
     24 #include <rsys/rsys.h>
     25 #include <rsys/str.h>
     26 #include <rsys/dynamic_array_int.h>
     27 #include <rsys/dynamic_array_double.h>
     28 #include <rsys/dynamic_array_char.h>
     29 #include <rsys/double3.h>
     30 #include <rsys/float3.h>
     31 
     32 #include <star/sg3d.h>
     33 #include <star/sg3d_sencXd_helper.h>
     34 #include <star/senc3d.h>
     35 
     36 #include <stdio.h>
     37 #include <stdlib.h>
     38 
     39 #define OKP(Expr) if((Expr) < 0) { res=RES_IO_ERR; goto error; }
     40 
     41 static int
     42 int_compare(const void *a_, const void *b_) {
     43   const int *a = (const int *)a_;
     44   const int *b = (const int *)b_;
     45   ASSERT(a && b);
     46   return (*a > *b) - (*a < *b);
     47 }
     48 
     49 static res_T
     50 write_ascii_stl
     51   (const char* filename,
     52    double* coords,
     53    const size_t count)
     54 {
     55   res_T res = RES_OK;
     56   size_t i;
     57   size_t actual_trg_count = 0;
     58   FILE* stl_file = NULL;
     59 
     60   ASSERT(filename && (coords || count == 0));
     61 
     62   stl_file = fopen(filename, "w");
     63   if(!stl_file) {
     64     res = RES_IO_ERR;
     65     goto error;
     66   }
     67 
     68   OKP(fprintf(stl_file, "solid %s\n", filename));
     69 
     70   for(i = 0; i < count; i++) {
     71     int k;
     72     double vtx[3][3], tmp[3], edge1[3], edge2[3], n[3] = { 0,0,0 };
     73     const double zero[3] = { 0,0,0 };
     74     d3_set(vtx[0], coords+9*i+0);
     75     d3_set(vtx[1], coords+9*i+3);
     76     d3_set(vtx[2], coords+9*i+6);
     77     d3_sub(edge1, vtx[1], vtx[0]);
     78     d3_sub(edge2, vtx[2], vtx[0]);
     79     d3_cross(tmp, edge1, edge2);
     80     if(d3_eq(tmp, zero))
     81       continue;
     82     d3_normalize(n, tmp);
     83 
     84     OKP(fprintf(stl_file, "  facet normal %g %g %g\n", SPLIT3(n)));
     85     OKP(fprintf(stl_file, "    outer loop\n"));
     86     for(k = 0; k < 3; k++) {
     87       OKP(fprintf(stl_file, "      vertex %.16g %.16g %.16g\n", SPLIT3(vtx[k])));
     88     }
     89     OKP(fprintf(stl_file, "    endloop\n"));
     90     OKP(fprintf(stl_file, "  endfacet\n"));
     91     actual_trg_count++;
     92   }
     93   OKP(fprintf(stl_file, "endsolid \n"));
     94   if(actual_trg_count != count) {
     95     size_t del_count = count - actual_trg_count;
     96     ASSERT(actual_trg_count < count);
     97     log_warning(get_device(),
     98         "In file '%s': %lu nul triangle(s) removed in exported geometry.\n",
     99         filename, (unsigned long)del_count);
    100   }
    101 
    102 exit:
    103     if(stl_file) fclose(stl_file);
    104   return res;
    105 error:
    106   log_error(get_device(), "Error: could not export to STL file '%s': %s\n",
    107     filename, res_to_cstr(res));
    108   goto exit;
    109 }
    110 
    111 static res_T
    112 write_binary_stl
    113   (const char* filename,
    114    double* coords,
    115    const size_t count)
    116 {
    117   res_T res = RES_OK;
    118   unsigned i;
    119   unsigned trg_count, actual_trg_count = 0;
    120   char header[80] = "Binary STL";
    121   FILE* stl_file = NULL;
    122   long int offset;
    123 
    124   ASSERT(filename && (coords || count == 0));
    125 
    126   stl_file = fopen(filename, "wb");
    127   if(!stl_file || count > UINT_MAX) {
    128     res = RES_IO_ERR;
    129     goto error;
    130   }
    131 
    132   if(1 != fwrite(header, sizeof(header), 1, stl_file)) {
    133     res = RES_IO_ERR;
    134     goto error;
    135   }
    136 
    137   trg_count = (unsigned)count;
    138   offset = ftell(stl_file);
    139   if(1 != fwrite(&trg_count, 4, 1, stl_file)) {
    140     res = RES_IO_ERR;
    141     goto error;
    142   }
    143 
    144   for(i = 0; i < trg_count; i++) {
    145     float tmp[3], edge1[3], edge2[3], zero[3] = { 0,0,0 };
    146     struct {
    147       float n[3];
    148       float vrtx[3][3];
    149       unsigned short attrib;
    150     } trg;
    151     f3_set_d3(trg.vrtx[0], coords+9*i+0);
    152     f3_set_d3(trg.vrtx[1], coords+9*i+3);
    153     f3_set_d3(trg.vrtx[2], coords+9*i+6);
    154     f3_sub(edge1, trg.vrtx[1], trg.vrtx[0]);
    155     f3_sub(edge2, trg.vrtx[2], trg.vrtx[0]);
    156     f3_cross(tmp, edge1, edge2);
    157     if(f3_eq(tmp, zero))
    158       continue;
    159     f3_normalize(trg.n, tmp);
    160     trg.attrib = 0;
    161     if(1 != fwrite(&trg, 50, 1, stl_file)) {
    162       res = RES_IO_ERR;
    163       goto error;
    164     }
    165     actual_trg_count++;
    166   }
    167   if(actual_trg_count != trg_count) {
    168     ASSERT(actual_trg_count < trg_count);
    169     log_warning(get_device(),
    170         "In file '%s': %u nul triangle(s) removed in exported geometry.\n",
    171         filename, trg_count - actual_trg_count);
    172     /* Need to change count */
    173     if(0 != fseek(stl_file, offset, SEEK_SET)) {
    174       res = RES_IO_ERR;
    175       goto error;
    176     }
    177     if(1 != fwrite(&actual_trg_count, 4, 1, stl_file)) {
    178       res = RES_IO_ERR;
    179       goto error;
    180     }
    181   }
    182 
    183 exit:
    184     if(stl_file) fclose(stl_file);
    185   return res;
    186 error:
    187   log_error(get_device(), "Error: could not export to STL file '%s': %s\n",
    188     filename, res_to_cstr(res));
    189   goto exit;
    190 }
    191 
    192 /* Accumulate the 2D tags of geometry (if it is 2D) or the tags of its boundary
    193  * (if it is 3D).
    194  * Trigger a call to sync_device() */
    195 res_T
    196 get_2d_tags
    197   (struct scad_geometry* geometry,
    198    struct darray_int* tags)
    199 {
    200   int* dimTags_to_free = NULL;
    201   size_t i;
    202   int ierr = 0;
    203   int* data;
    204   size_t sz;
    205   res_T res = RES_OK;
    206 
    207   ASSERT(geometry && tags);
    208 
    209   sz = geometry->gmsh_dimTags_n;
    210   data = geometry->gmsh_dimTags;
    211   ASSERT(sz % 2 == 0);
    212 
    213   ERR(sync_device());
    214   for(i = 0; i < sz; i += 2) {
    215     int dim = data[i];
    216     int tag = data[i+1];
    217     size_t j;
    218     if(dim == 3) {
    219       int* face_dimTags = NULL;
    220       size_t face_dimTags_n, count;
    221       /* Get 2d components of the 3d object */
    222       gmshModelGetBoundary(data+i, 2, &face_dimTags, &face_dimTags_n, 1, 0, 0, &ierr);
    223       dimTags_to_free = face_dimTags;
    224       ERR(gmsh_err_to_res_T(ierr));
    225       ASSERT(face_dimTags_n % 2 == 0);
    226       /* copy tags in the output */
    227       count = darray_int_size_get(tags);
    228       ERR(darray_int_reserve(tags, count + face_dimTags_n/2));
    229       for(j = 0; j < face_dimTags_n; j += 2) {
    230         ASSERT(face_dimTags[j] == 2); /* dim == 2 */
    231         ERR(darray_int_push_back(tags, &face_dimTags[j+1]));
    232       }
    233       gmshFree(dimTags_to_free); dimTags_to_free = NULL;
    234     }
    235     else if(dim == 2) {
    236       ERR(darray_int_push_back(tags, &tag));
    237     } else {
    238       res = RES_BAD_ARG;
    239       goto error;
    240     }
    241   }
    242 exit:
    243   gmshFree(dimTags_to_free);
    244   return res;
    245 error:
    246   goto exit;
    247 }
    248 
    249 /******************************************************************************
    250  * Exported functions
    251  *****************************************************************************/
    252 res_T
    253 scad_synchronize
    254   (void)
    255 {
    256   res_T res = RES_OK;
    257   struct scad_device* dev = get_device();
    258   int ierr;
    259 
    260   /* Cannot use check_device to avoid an infinite loop! */
    261   if(!dev) {
    262     /* No logger available for a message */
    263     fprintf(stderr,
    264        "%s: cannot call API functions if star-cad is not initialized.\n",
    265        FUNC_NAME);
    266     res = RES_BAD_ARG;
    267     goto error;
    268   }
    269 
    270   gmshModelOccSynchronize(&ierr);
    271   dev->need_synchro = dev->options.Misc.DebugAutoSync;
    272   ERR(gmsh_err_to_res_T(ierr));
    273 
    274 exit:
    275   return res;
    276 error:
    277   goto exit;
    278 }
    279 
    280 res_T
    281 scad_run_ui(void)
    282 {
    283   res_T res = RES_OK;
    284   struct scad_device* dev = get_device();
    285   int ierr;
    286 
    287   /* Cannot use check_device to avoid an infinite loop! */
    288   if(!dev) {
    289     /* No logger available for a message */
    290     fprintf(stderr,
    291        "%s: cannot call API functions if star-cad is not initialized.\n",
    292        FUNC_NAME);
    293     res = RES_BAD_ARG;
    294     goto error;
    295   }
    296 
    297   if(dev->options.Misc.SynchronizeOnRunUI) {
    298     ERR(sync_device());
    299   }
    300 
    301   gmshFltkRun(&ierr);
    302   if(ierr) {
    303     log_error(dev, "Cannot call FLTK: you probably need to build gmsh locally.\n");
    304   }
    305   ERR(gmsh_err_to_res_T(ierr));
    306 
    307 exit:
    308   return res;
    309 error:
    310   goto exit;
    311 }
    312 
    313 size_t
    314 scad_get_dimtag_refcount
    315   (const int dim,
    316    const int tag)
    317 {
    318   struct tag_desc* desc = NULL;
    319 
    320   if(check_device(FUNC_NAME) != RES_OK) return SIZE_MAX;
    321   desc = device_get_description(dim, tag);
    322   if(!desc) return SIZE_MAX;
    323 
    324   return desc->refcount;
    325 }
    326 
    327 res_T
    328 scad_scene_write
    329   (const char* name)
    330 {
    331   int ierr;
    332   res_T res = RES_OK;
    333 
    334   if(!name) {
    335     res = RES_BAD_ARG;
    336     goto error;
    337   }
    338 
    339   ERR(check_device(FUNC_NAME));
    340   ERR(sync_device());
    341   gmshWrite(name, &ierr);
    342   ERR(gmsh_err_to_res_T(ierr));
    343 
    344 exit:
    345   return res;
    346 error:
    347   goto exit;
    348 }
    349 
    350 res_T
    351 scad_stl_get_data
    352   (struct scad_geometry* geometry,
    353    struct darray_double* triangles)
    354 {
    355   return scad_stl_get_data_partial(geometry, NULL, 0, triangles);
    356 }
    357 
    358 res_T
    359 scad_stl_get_data_partial
    360   (struct scad_geometry* geometry,
    361    struct scad_geometry** dont,
    362    const size_t dont_count,
    363    struct darray_double* triangles)
    364 {
    365   size_t* nodeTags = NULL;
    366   double* coord = NULL;
    367   double* pCoord = NULL;
    368   size_t i, count0, tcount, sz, dsz = 0;
    369   int ierr = 0;
    370   int *data, *ddata = NULL;
    371   struct scad_device* dev = get_device();
    372   struct mem_allocator* allocator = NULL;
    373   struct darray_int tags;
    374   struct darray_int dtags;
    375   int tags_initialized = 0;
    376   res_T res = RES_OK;
    377 
    378   if(!geometry || !triangles || (!dont && dont_count > 0)) {
    379     res = RES_BAD_ARG;
    380     goto error;
    381   }
    382 
    383   ERR(check_device(FUNC_NAME));
    384 
    385   allocator = dev->allocator;
    386   darray_int_init(allocator, &tags);
    387   darray_int_init(allocator, &dtags);
    388   tags_initialized = 1;
    389 
    390   /* Get 2d tags */
    391   ERR(get_2d_tags(geometry, &tags));
    392   sz = darray_int_size_get(&tags);
    393   data = darray_int_data_get(&tags);
    394 
    395   if(dont) {
    396     for(i = 0; i < dont_count; i++) {
    397       ERR(get_2d_tags(dont[i], &dtags));
    398     }
    399     dsz = darray_int_size_get(&dtags);
    400     ddata = darray_int_data_get(&dtags);
    401     /* Sort ddata tags to enable bsearch */
    402     qsort(ddata, dsz, sizeof(*ddata), int_compare);
    403   }
    404 
    405   /* Collect triangles */
    406   count0 = tcount = darray_double_size_get(triangles);
    407   for(i = 0; i < sz; i++) {
    408     size_t j, coord_n, pCoord_n, nodeTags_n;
    409     const int type = 2; /* 3-node triangles (see src/common/GmshDefines.h) */
    410     const int tag = data[i];
    411     if(bsearch(&tag, ddata, dsz, sizeof(*ddata), int_compare)) {
    412       /* this tag is part of the dont list: don't collect */
    413       continue;
    414     }
    415     gmshModelMeshGetNodesByElementType(type, &nodeTags, &nodeTags_n,
    416         &coord, &coord_n, &pCoord, &pCoord_n, tag, 0, &ierr);
    417     ERR(gmsh_err_to_res_T(ierr));
    418     tcount += coord_n;
    419     ERR(darray_double_reserve(triangles, tcount));
    420     for(j = 0; j < coord_n; j += 9) {
    421 #define PUSH3(A, D) { \
    422   ERR(darray_double_push_back((A), (D)+0)); \
    423   ERR(darray_double_push_back((A), (D)+1)); \
    424   ERR(darray_double_push_back((A), (D)+2)); \
    425 }
    426       /* Keep the triangle; if its normal is reversed, change vertices' order */
    427       PUSH3(triangles, coord+j+0);
    428       PUSH3(triangles, coord+j+3);
    429       PUSH3(triangles, coord+j+6);
    430     }
    431 #undef PUSH3
    432     gmshFree(nodeTags); nodeTags = NULL;
    433     gmshFree(coord); coord = NULL;
    434     gmshFree(pCoord); pCoord = NULL;
    435   }
    436   ASSERT(tcount == darray_double_size_get(triangles));
    437   if(count0 == tcount) {
    438     if(str_is_empty(&geometry->name)) {
    439       log_message(dev, "No triangle collected for unnamed geometry '%p'.\n",
    440           (void*)geometry);
    441     } else {
    442       log_message(dev, "No triangle collected for geometry '%s'.\n",
    443           str_cget(&geometry->name));
    444     }
    445   }
    446 
    447 exit:
    448   if(!allocator) {
    449     /* Early error, nothing was allocated */
    450     return res;
    451   }
    452   if(tags_initialized) {
    453     darray_int_release(&tags);
    454     darray_int_release(&dtags);
    455   }
    456   gmshFree(nodeTags);
    457   gmshFree(coord);
    458   gmshFree(pCoord);
    459   return res;
    460 error:
    461   if(dev) {
    462     log_error(dev, "%s: could not get STL data -- %s\n",
    463       FUNC_NAME, res_to_cstr(res));
    464   }
    465   goto exit;
    466 }
    467 
    468 static INLINE void
    469 get_indices(const unsigned itri, unsigned ids[3], void* context)
    470 {
    471   const double* ctx = context;
    472   ASSERT(ids && ctx); (void)ctx;
    473   ids[0] = itri * 3 + 0;
    474   ids[1] = itri * 3 + 1;
    475   ids[2] = itri * 3 + 2;
    476 }
    477 
    478 static INLINE void
    479 get_position(const unsigned ivert, double pos[3], void* context)
    480 {
    481   const double* ctx = context;
    482   ASSERT(pos && ctx);
    483   pos[0] = ctx[ivert * 3 + 0];
    484   pos[1] = ctx[ivert * 3 + 1];
    485   pos[2] = ctx[ivert * 3 + 2];
    486 }
    487 
    488 static res_T
    489 scad_stl_sort_orientation
    490   (struct darray_double* triangles,
    491    const char* set_name,
    492    const enum scad_normals_orientation orientation)
    493 {
    494   res_T res = RES_OK;
    495   struct scad_device* dev = get_device();
    496   struct mem_allocator* allocator = dev->allocator;
    497   struct logger* logger = dev->logger;
    498   double* coord;
    499   size_t coord_n;
    500   unsigned i, ecount, tcount_in, tcount_out, vcount_in, utcount_in, vcount;
    501   struct sg3d_device* sg3d = NULL;
    502   struct sg3d_geometry* geom = NULL;
    503   struct senc3d_device* senc3d = NULL;
    504   struct senc3d_scene* senc3d_scn = NULL;
    505   struct senc3d_enclosure* enclosure = NULL;
    506   struct senc3d_enclosure_header header;
    507   int convention, initialized = 0;
    508   struct sg3d_geometry_add_callbacks callbacks = SG3D_ADD_CALLBACKS_NULL__;
    509   struct darray_double new_triangles;
    510   char* tin = NULL;
    511   char *contact_0 = NULL;
    512   unsigned ocount = 0;
    513   struct darray_double dblsided;
    514   struct str dbl_name;
    515   int dblsided_initialized = 0;
    516 
    517   if(!triangles || !set_name) {
    518     res = RES_BAD_ARG;
    519     goto error;
    520   }
    521 
    522   coord_n = darray_double_size_get(triangles);
    523   if(coord_n % 9 || coord_n > UINT_MAX) {
    524     res = RES_BAD_ARG;
    525     goto error;
    526   }
    527   if(coord_n == 0 || orientation == SCAD_KEEP_NORMALS_UNCHANGED) {
    528     goto exit;
    529   }
    530 
    531   tcount_in = (unsigned)(coord_n/9);
    532   vcount_in = tcount_in * 3;
    533   coord = darray_double_data_get(triangles);
    534 
    535   ERR(sg3d_device_create(logger, allocator, dev->verbose, &sg3d));
    536   ERR(sg3d_geometry_create(sg3d, &geom));
    537 
    538   callbacks.get_indices = get_indices;
    539   callbacks.get_position = get_position;
    540 
    541   ERR(sg3d_geometry_reserve(geom, vcount_in, tcount_in, 0));
    542   ERR(sg3d_geometry_add(geom, vcount_in, tcount_in, &callbacks, coord));
    543 
    544   ERR(sg3d_geometry_get_unique_vertices_count(geom, &vcount));
    545   ERR(sg3d_geometry_get_unique_triangles_count(geom, &utcount_in));
    546   if(utcount_in != tcount_in) {
    547     ASSERT(tcount_in > utcount_in);
    548     log_warning(get_device(),
    549         "Triangles duplicates found when sorting out normals (%u / %u) in set '%s'.\n",
    550         tcount_in - utcount_in, tcount_in, set_name);
    551   }
    552   if(orientation == SCAD_FORCE_NORMALS_OUTWARD)
    553    convention = SENC3D_CONVENTION_NORMAL_BACK | SENC3D_CONVENTION_NORMAL_OUTSIDE;
    554   else
    555    convention = SENC3D_CONVENTION_NORMAL_BACK | SENC3D_CONVENTION_NORMAL_INSIDE;
    556 
    557   ERR(senc3d_device_create(logger, allocator,
    558         dev->options.Geometry.OCCParallel ? SENC3D_NTHREADS_DEFAULT : 1,
    559         dev->verbose, &senc3d));
    560   ERR(senc3d_scene_create(senc3d, convention,
    561         utcount_in, sg3d_sencXd_geometry_get_indices, NULL,
    562         vcount, sg3d_sencXd_geometry_get_position, geom, &senc3d_scn));
    563 
    564   ERR(senc3d_scene_get_overlapping_triangles_count(senc3d_scn, &ocount));
    565   if(ocount) {
    566     logger_print(logger, LOG_ERROR,
    567       "Overlapping triangles found when sorting out normals (%u / %u) "
    568         "in set '%s'.\n",
    569     utcount_in - ocount, utcount_in, set_name);
    570     res = RES_BAD_ARG;
    571     goto error;
    572   }
    573 
    574   ERR(senc3d_scene_get_enclosure_count(senc3d_scn, &ecount));
    575   /* If there is at least 1 triangle, there is at least enclosure #0 that
    576    * represents the outside of the whole scene.
    577    * For the model to allow forcing normals, we need a model defining an
    578    * inside/outside frontier. These models have at least 2 enclosures. */
    579   if(ecount < 2) {
    580     /* Must define a closed volume */
    581     log_error(get_device(),
    582         "Triangle set '%s' doesn't define a closed volume.\n", set_name);
    583     res = RES_BAD_ARG;
    584     goto error;
    585   } else {
    586     /* For the model to allow forcing normals, we reject models defining
    587      * enclosures included inside other enclosures and models with triangles
    588      * that have their 2 sides in the output, as their normal's orientation
    589      * could not satisfy the requirement of the 2 sides.
    590      * To properly detect enclosures inside enclosures, we rely on the number of
    591      * output triangles, that must be identical to the number of input
    592      * triangles. */
    593     unsigned e, enclosures[2];
    594     size_t two = 0;
    595     darray_double_init(allocator, &new_triangles);
    596     initialized = 1;
    597     contact_0 = MEM_CALLOC(allocator, ecount, sizeof(*contact_0));
    598     if(!contact_0) {
    599       res = RES_MEM_ERR;
    600       goto error;
    601     }
    602     /* Triangles of enclosure #0 are also part of other enclosures. To have them
    603      * with the expected normal orientation we have to get them from these other
    604      * enclosures. */
    605     for(e = 1; e < ecount; e++) {
    606       ERR(senc3d_scene_get_enclosure(senc3d_scn, e, &enclosure));
    607       ERR(senc3d_enclosure_get_header(enclosure, &header));
    608       two += (header.primitives_count - header.unique_primitives_count);
    609       for(i = 0; i < header.unique_primitives_count; i++) {
    610         enum senc3d_side side;
    611         unsigned gid;
    612         ERR(senc3d_enclosure_get_triangle_id(enclosure, i, &gid, &side));
    613         ERR(senc3d_scene_get_triangle_enclosures(senc3d_scn, gid, enclosures));
    614         if(enclosures[0] == 0 || enclosures[1] == 0) {
    615           contact_0[e] = 1;
    616           break;
    617         }
    618       }
    619       ERR(senc3d_enclosure_ref_put(enclosure));
    620       enclosure = NULL;
    621     }
    622     if(two > 0) {
    623       size_t idx;
    624       res_T r;
    625       log_error(get_device(),
    626           "Triangle set '%s' defines an invalid closed volume "
    627           "(%u / %u triangles with both sides in).\n",
    628           set_name, (unsigned)two, utcount_in);
    629       /* Output the two-sides-in triangles */
    630       darray_double_init(allocator, &dblsided);
    631       str_init(allocator, &dbl_name);
    632       dblsided_initialized = 1;
    633       ERR(darray_double_reserve(&dblsided, 9 * two));
    634       for(e = 0; e < ecount; e++) {
    635         ERR(senc3d_scene_get_enclosure(senc3d_scn, e, &enclosure));
    636         ERR(senc3d_enclosure_get_header(enclosure, &header));
    637         for(i = header.unique_primitives_count; i < header.primitives_count; i++) {
    638           enum senc3d_side side;
    639           unsigned gid, j, k;
    640           unsigned trg[3];
    641           double v[3];
    642           ERR(senc3d_enclosure_get_triangle_id(enclosure, i, &gid, &side));
    643           ERR(senc3d_enclosure_get_triangle(enclosure, gid, trg));
    644           ASSERT(side == (SENC3D_FRONT | SENC3D_BACK));
    645           for(j = 0; j < 3; j++) {
    646             ERR(senc3d_enclosure_get_vertex(enclosure, trg[j], v));
    647             for(k = 0; k < 3; k++) {
    648               ERR(darray_double_push_back(&dblsided, v+k));
    649             }
    650           }
    651         }
    652         ERR(senc3d_enclosure_ref_put(enclosure));
    653         enclosure = NULL;
    654       }
    655       idx = strlen(set_name);
    656       ASSERT(idx > 3 && 0 == strcmp(set_name + idx - 4, ".stl"));
    657       str_set(&dbl_name, set_name);
    658       str_insert(&dbl_name, idx - 4, "_double_sided_triangles");
    659       r = scad_stl_data_write(&dblsided, str_cget(&dbl_name),
    660           SCAD_KEEP_NORMALS_UNCHANGED, 0);
    661       if(r == RES_OK) {
    662         log_error(get_device(),
    663             "Saved double sided triangles to file '%s'.\n",
    664             str_cget(&dbl_name));
    665       } else {
    666         log_error(get_device(),
    667             "Could not save double sided triangles to a file.\n");
    668       }
    669       res = RES_BAD_ARG;
    670       goto error;
    671     }
    672     ERR(darray_double_reserve(&new_triangles, 9 * utcount_in));
    673     tin = MEM_CALLOC(allocator, 1, 9 * utcount_in);
    674     if(!tin) {
    675       res= RES_MEM_ERR;
    676       goto error;
    677     }
    678     for(e = 1; e < ecount; e++) {
    679       if(!contact_0[e]) {
    680         /* Output only enclosures in contact with enclosure #0 */
    681         continue;
    682       }
    683       ERR(senc3d_scene_get_enclosure(senc3d_scn, e, &enclosure));
    684       ERR(senc3d_enclosure_get_header(enclosure, &header));
    685       for(i = 0; i < header.unique_primitives_count; i++) {
    686         unsigned n, k, idx, vrt[3];
    687         enum senc3d_side side;
    688         /* Ensure that input triangles are included only once */
    689         ERR(senc3d_enclosure_get_triangle_id(enclosure, i, &idx, &side));
    690         if(tin[idx]) {
    691           /* Allready in output */
    692           char s = (side == SENC3D_FRONT) ? 1 : 2;;
    693           if(s == tin[idx])
    694             continue; /* Same side => OK (should not happen?) */
    695           log_error(get_device(),
    696               "Triangle set '%s' defines a topology that doesn't allow forcing normals"
    697               " (found triangle(s) with 2 sides in).\n",
    698               set_name);
    699           res = RES_BAD_ARG;
    700           goto error;
    701         }
    702         tin[idx] = (side == SENC3D_FRONT) ? 1 : 2;;
    703         /* Get the vertices as they are ordered in the enclosure's mesh */
    704         ERR(senc3d_enclosure_get_triangle(enclosure, i, vrt));
    705         /* Rewrite vertices according to enclosure's mesh */
    706         for(n = 0; n < 3; n++) {
    707           double pos[3];
    708           ERR(senc3d_enclosure_get_vertex(enclosure, vrt[n], pos));
    709           for(k = 0; k < 3; k++) {
    710             ERR(darray_double_push_back(&new_triangles, pos+k));
    711           }
    712         }
    713       }
    714       ERR(senc3d_enclosure_ref_put(enclosure));
    715       enclosure = NULL;
    716     }
    717     tcount_out = (unsigned)(darray_double_size_get(&new_triangles)/9);
    718     ASSERT(utcount_in >= tcount_out);
    719     if(utcount_in != tcount_out) {
    720       log_error(get_device(),
    721           "Triangle set '%s' defines a topology that doesn't allow forcing normals"
    722           " (%u /%u triangles not part of the output).\n",
    723           set_name, utcount_in - tcount_out, utcount_in);
    724       res = RES_BAD_ARG;
    725       goto error;
    726     }
    727     ERR(darray_double_copy_and_release(triangles, &new_triangles));
    728     initialized = 0;
    729   }
    730 
    731 exit:
    732   MEM_RM(allocator, contact_0);
    733   MEM_RM(allocator, tin);
    734   if(dblsided_initialized) {
    735     darray_double_release(&dblsided);
    736     str_release(&dbl_name);
    737   }
    738   if(initialized) darray_double_release(&new_triangles);
    739   if(sg3d) SG3D(device_ref_put(sg3d));
    740   if(geom) SG3D(geometry_ref_put(geom));
    741   if(senc3d) SENC3D(device_ref_put(senc3d));
    742   if(senc3d_scn) SENC3D(scene_ref_put(senc3d_scn));
    743   if(enclosure) SENC3D(enclosure_ref_put(enclosure));
    744   return res;
    745 error:
    746   goto exit;
    747 }
    748 
    749 res_T
    750 scad_stl_data_write
    751   (const struct darray_double* triangles,
    752    const char* filename,
    753    const enum scad_normals_orientation orientation,
    754    const int binary)
    755 {
    756   res_T res = RES_OK;
    757   res_T tmp_res = RES_OK;
    758   struct scad_device* dev = get_device();
    759   struct darray_double sorted;
    760   int initialized = 0;
    761   double* coord;
    762   size_t coord_n;
    763 
    764   if(!triangles || !filename) {
    765     res = RES_BAD_ARG;
    766     goto error;
    767   }
    768 
    769   ERR(check_device(FUNC_NAME));
    770 
    771   coord_n = darray_double_size_get(triangles);
    772   if(coord_n % 9) {
    773     res = RES_BAD_ARG;
    774     goto error;
    775   }
    776 
    777   darray_double_init(dev->allocator, &sorted);
    778   initialized = 1;
    779   ERR(darray_double_copy(&sorted, triangles));
    780 
    781   /* If sort_orientation fails, try to write the file anyway to allow debugging */
    782   tmp_res = scad_stl_sort_orientation(&sorted, filename, orientation);
    783   coord = darray_double_data_get(&sorted);
    784   coord_n = darray_double_size_get(&sorted);
    785   if(binary) ERR(write_binary_stl(filename, coord, coord_n/9));
    786   else ERR(write_ascii_stl(filename, coord, coord_n/9));
    787   ERR(tmp_res);
    788 
    789 exit:
    790   if(initialized) darray_double_release(&sorted);
    791   return res;
    792 error:
    793   goto exit;
    794 }
    795 
    796 res_T
    797 scad_stl_export
    798   (struct scad_geometry* geometry,
    799    const char* file_name,
    800    const enum scad_normals_orientation orientation,
    801    const int binary)
    802 {
    803   res_T res = RES_OK;
    804   struct darray_double trg;
    805   struct scad_device* dev = get_device();
    806   struct str filename;
    807   int initialized = 0;
    808 
    809   if(!geometry || (!file_name && str_is_empty(&geometry->name)) ||
    810       (orientation != SCAD_KEEP_NORMALS_UNCHANGED &&
    811        orientation != SCAD_FORCE_NORMALS_OUTWARD &&
    812        orientation != SCAD_FORCE_NORMALS_INWARD))
    813   {
    814     res = RES_BAD_ARG;
    815     goto error;
    816   }
    817 
    818   ERR(check_device(FUNC_NAME));
    819 
    820   darray_double_init(dev->allocator, &trg);
    821   str_init(dev->allocator, &filename);
    822   initialized = 1;
    823   if(file_name) {
    824     ERR(str_set(&filename, file_name));
    825   } else {
    826     if(str_is_empty(&geometry->name)) {
    827       res = RES_BAD_ARG;
    828       goto error;
    829     }
    830     ERR(str_copy(&filename, &geometry->name));
    831   }
    832   ERR(str_append(&filename, ".stl"));
    833   ERR(scad_stl_get_data(geometry, &trg));
    834   ERR(scad_stl_data_write(&trg, str_cget(&filename), orientation, binary));
    835 
    836 exit:
    837   if(initialized) {
    838     darray_double_release(&trg);
    839     str_release(&filename);
    840   }
    841   return res;
    842 error:
    843   goto exit;
    844 }
    845 
    846 res_T
    847 scad_stl_export_partial
    848   (struct scad_geometry* geometry,
    849    struct scad_geometry** dont,
    850    const size_t dont_count,
    851    const char* file_name,
    852    const enum scad_normals_orientation orientation,
    853    const int binary)
    854 {
    855   res_T res = RES_OK;
    856   struct darray_double trg;
    857   struct scad_device* dev = get_device();
    858   struct str filename;
    859   int initialized = 0;
    860 
    861   if(!geometry || (!file_name && str_is_empty(&geometry->name)) ||
    862       (!dont && dont_count > 0) ||
    863       (orientation != SCAD_KEEP_NORMALS_UNCHANGED &&
    864        orientation != SCAD_FORCE_NORMALS_OUTWARD &&
    865        orientation != SCAD_FORCE_NORMALS_INWARD))
    866   {
    867     res = RES_BAD_ARG;
    868     goto error;
    869   }
    870 
    871   ERR(check_device(FUNC_NAME));
    872 
    873   darray_double_init(dev->allocator, &trg);
    874   str_init(dev->allocator, &filename);
    875   initialized = 1;
    876   if(file_name) {
    877     ERR(str_set(&filename, file_name));
    878   } else {
    879     if(str_is_empty(&geometry->name)) {
    880       res = RES_BAD_ARG;
    881       goto error;
    882     }
    883     ERR(str_copy(&filename, &geometry->name));
    884   }
    885   ERR(str_append(&filename, ".stl"));
    886   ERR(scad_stl_get_data_partial(geometry, dont, dont_count, &trg));
    887   ERR(scad_stl_data_write(&trg, str_cget(&filename), orientation, binary));
    888 
    889 exit:
    890   if(initialized) {
    891     darray_double_release(&trg);
    892     str_release(&filename);
    893   }
    894   return res;
    895 error:
    896   goto exit;
    897 }
    898 
    899 res_T
    900 scad_stl_export_split
    901   (struct scad_geometry* geometry,
    902    const char* file_name,
    903    const enum scad_normals_orientation orientation,
    904    const int binary)
    905 {
    906   size_t i;
    907   struct scad_geometry** parts = NULL;
    908   struct scad_device* dev = get_device();
    909   size_t cpt = 0, count;
    910   struct str name;
    911   int init = 0;
    912   res_T res = RES_OK;
    913 
    914   if(!geometry || (!file_name && str_is_empty(&geometry->name))) {
    915     res = RES_BAD_ARG;
    916     goto error;
    917   }
    918 
    919   ERR(check_device(FUNC_NAME));
    920 
    921   ERR(scad_geometry_explode(geometry, &parts, &count));
    922   ASSERT(count*2 == geometry->gmsh_dimTags_n);
    923   str_init(dev->allocator, &name);
    924   init = 1;
    925   for(i = 0; i < count; i++) {
    926     ERR(str_printf(&name, "%s_%ld", file_name, cpt++));
    927     ERR(scad_stl_export(parts[i], str_cget(&name), orientation, binary));
    928   }
    929 
    930 exit:
    931   if(parts) {
    932     for(i = 0; i < count; i++) SCAD(geometry_ref_put(parts[i]));
    933     MEM_RM(dev->allocator, parts);
    934   }
    935   if(init) str_release(&name);
    936   return res;
    937 error:
    938   goto exit;
    939 }
    940 
    941 res_T
    942 scad_scene_mesh
    943   (void)
    944 {
    945   int ierr = 0;
    946   res_T res = RES_OK;
    947 
    948   ERR(check_device(FUNC_NAME));
    949 
    950   ERR(sync_device());
    951   gmshModelMeshGenerate(2, &ierr);
    952   ERR(gmsh_err_to_res_T(ierr));
    953 
    954 exit:
    955   return res;
    956 error:
    957   goto exit;
    958 }