git.mutantstargoat.com Git - meshfrac/blob - src/frac.c

   1 #include <stdlib.h>
   2 #include <string.h>
   3 #include <assert.h>
   4 #include "frac.h"
   5 #include "dynarr.h"
   6
   7 static void destroy_cell(struct frac_cell *cell);
   8 static int build_cell(struct fracture *frac, int cellidx);
   9
  10 int frac_init(struct fracture *frac)
  11 {
  12         frac->mesh = 0;
  13
  14         if(!(frac->cells = dynarr_alloc(0, sizeof *frac->cells))) {
  15                 return -1;
  16         }
  17         return 0;
  18 }
  19
  20 void frac_destroy(struct fracture *frac)
  21 {
  22         int i, num;
  23
  24         if(!frac) return;
  25
  26         if(frac->cells) {
  27                 num = dynarr_size(frac->cells);
  28                 for(i=0; i<num; i++) {
  29                         destroy_cell(frac->cells + i);
  30                 }
  31                 dynarr_free(frac->cells);
  32         }
  33 }
  34
  35 static int init_cell(struct frac_cell *cell)
  36 {
  37         memset(cell, 0, sizeof *cell);
  38         if(!(cell->mesh = cmesh_alloc())) {
  39                 return -1;
  40         }
  41         return 0;
  42 }
  43
  44 static void destroy_cell(struct frac_cell *cell)
  45 {
  46         if(!cell) return;
  47         free(cell->polys);
  48         cmesh_free(cell->mesh);
  49 }
  50
  51 void frac_mesh(struct fracture *frac, const struct cmesh *m)
  52 {
  53         frac->mesh = (struct cmesh*)m;
  54 }
  55
  56 int frac_point(struct fracture *frac, float x, float y, float z)
  57 {
  58         struct frac_cell cell;
  59         struct frac_cell *tmp;
  60
  61         init_cell(&cell);
  62         cgm_vcons(&cell.pt, x, y, z);
  63         if(!(tmp = dynarr_push(frac->cells, &cell))) {
  64                 destroy_cell(&cell);
  65                 return -1;
  66         }
  67         frac->cells = tmp;
  68         return 0;
  69 }
  70
  71 int frac_num_cells(struct fracture *frac)
  72 {
  73         return dynarr_size(frac->cells);
  74 }
  75
  76 /* --- step 1: generate a bunch of points (or let the user add them manually) */
  77
  78 int frac_gen_points(struct fracture *frac, int num)
  79 {
  80         int i;
  81         float x, y, z;
  82         cgm_vec3 bbmin, bbmax, delta;
  83
  84         if(!frac || !frac->mesh) return -1;
  85         if(!cmesh_poly_count(frac->mesh)) {
  86                 return -1;
  87         }
  88
  89         cmesh_aabbox(frac->mesh, &bbmin, &bbmax);
  90         delta = bbmax;
  91         cgm_vsub(&delta, &bbmin);
  92
  93         for(i=0; i<num; i++) {
  94                 x = (float)rand() / RAND_MAX * delta.x + bbmin.x;
  95                 y = (float)rand() / RAND_MAX * delta.y + bbmin.y;
  96                 z = (float)rand() / RAND_MAX * delta.z + bbmin.z;
  97
  98                 if(frac_point(frac, x, y, z) == -1) {
  99                         return -1;
 100                 }
 101         }
 102         return 0;
 103 }
 104
 105
 106 /* --- step 2: construct voronoi cells bounded by planes */
 107
 108 int frac_build_cells(struct fracture *frac)
 109 {
 110         int i;
 111
 112         for(i=0; i<dynarr_size(frac->cells); i++) {
 113                 if(build_cell(frac, i) == -1) {
 114                         return -1;
 115                 }
 116         }
 117
 118         return 0;
 119 }
 120
 121 static int build_cell(struct fracture *frac, int cellidx)
 122 {
 123         int i, j, num, clipres;
 124         int *valid_planes;
 125         struct plane plane;
 126         struct frac_cell *cell = frac->cells + cellidx;
 127         struct poly poly, clipped, *polys, *pptr;
 128         float bsize;
 129
 130         if(!(polys = dynarr_alloc(0, sizeof *polys))) {
 131                 return -1;
 132         }
 133
 134         cmesh_bsphere(frac->mesh, 0, &bsize);
 135         bsize *= 8;
 136
 137         num = dynarr_size(frac->cells);
 138         for(i=0; i<num; i++) {
 139                 if(i == cellidx) continue;
 140                 midplane(&plane, &frac->cells[i].pt, &cell->pt, frac->cell_gap);
 141                 if(plane_sdist(&plane, &cell->pt) > 0.0f) {
 142                         plane_poly(&plane, &poly, bsize);
 143                         if(!(pptr = dynarr_push(polys, &poly))) {
 144                                 return -1;
 145                         }
 146                         polys = pptr;
 147                 }
 148         }
 149
 150         num = dynarr_size(polys);
 151         valid_planes = alloca(num * sizeof *valid_planes);
 152         memset(valid_planes, 0xff, num * sizeof *valid_planes);
 153
 154         /* clip all planes against each other to end up with a convex cell */
 155         cell->num_polys = num;
 156         for(i=0; i<num; i++) {
 157                 for(j=0; j<num; j++) {
 158                         if(i == j || !valid_planes[j]) {
 159                                 continue;
 160                         }
 161
 162                         /* clip plane polygon i with plane j */
 163                         poly_plane(polys + j, &plane);
 164                         init_poly(&clipped);
 165                         if((clipres = clip_poly(&clipped, polys + i, &plane)) < 0) {
 166                                 /* completely clipped, mark it for removal */
 167                                 valid_planes[i] = 0;
 168                                 cell->num_polys--;
 169                                 destroy_poly(&clipped);
 170                                 break;
 171                         }
 172                         destroy_poly(polys + i);
 173                         polys[i] = clipped;
 174                 }
 175         }
 176
 177         if(!(cell->polys = malloc(cell->num_polys * sizeof *cell->polys))) {
 178                 return -1;
 179         }
 180         pptr = cell->polys;
 181         for(i=0; i<num; i++) {
 182                 if(valid_planes[i]) {
 183                         assert(pptr - cell->polys < cell->num_polys);
 184                         *pptr++ = polys[i];
 185                 } else {
 186                         destroy_poly(polys + i);
 187                 }
 188         }
 189         dynarr_free(polys);
 190         return 0;
 191 }
 192
 193 static int mesh_poly(struct poly *poly, const struct cmesh *mesh, int faceidx)
 194 {
 195         int i, vsz, nsz, uvsz;
 196         struct vertex *tmpvert, vert = {0};
 197         const float *varr, *narr, *uvarr;
 198         unsigned int vidx;
 199
 200         if(init_poly(poly) == -1) {
 201                 return -1;
 202         }
 203
 204         varr = cmesh_attrib_ro(mesh, CMESH_ATTR_VERTEX);
 205         narr = cmesh_attrib_ro(mesh, CMESH_ATTR_NORMAL);
 206         uvarr = cmesh_attrib_ro(mesh, CMESH_ATTR_TEXCOORD);
 207         vsz = cmesh_attrib_nelem(mesh, CMESH_ATTR_VERTEX);
 208         nsz = cmesh_attrib_nelem(mesh, CMESH_ATTR_NORMAL);
 209         uvsz = cmesh_attrib_nelem(mesh, CMESH_ATTR_TEXCOORD);
 210
 211         for(i=0; i<3; i++) {
 212                 if(cmesh_indexed(mesh)) {
 213                         vidx = cmesh_index_ro(mesh)[faceidx * 3 + i];
 214                 } else {
 215                         vidx = faceidx * 3 + i;
 216                 }
 217                 vert.pos.x = varr[vidx * vsz];
 218                 vert.pos.y = varr[vidx * vsz + 1];
 219                 vert.pos.z = varr[vidx * vsz + 2];
 220                 if(narr) {
 221                         vert.norm.x = narr[vidx * nsz];
 222                         vert.norm.y = narr[vidx * nsz + 1];
 223                         vert.norm.z = narr[vidx * nsz + 2];
 224                 }
 225                 if(uvarr) {
 226                         vert.uv.x = uvarr[vidx * uvsz];
 227                         vert.uv.y = uvarr[vidx * uvsz + 1];
 228                 }
 229
 230                 if(!(tmpvert = dynarr_push(poly->verts, &vert))) {
 231                         destroy_poly(poly);
 232                         return -1;
 233                 }
 234                 poly->verts = tmpvert;
 235         }
 236         return 0;
 237 }
 238
 239 #define ADD_VERTEX(mesh, vert) \
 240         do { \
 241                 cmesh_normal(mesh, (vert)->norm.x, (vert)->norm.y, (vert)->norm.z); \
 242                 cmesh_texcoord(mesh, (vert)->uv.x, (vert)->uv.y, 0); \
 243                 if(cmesh_vertex(mesh, (vert)->pos.x, (vert)->pos.y, (vert)->pos.z) == -1) { \
 244                         fprintf(stderr, "SKATA\n"); \
 245                         abort(); \
 246                 } \
 247         } while(0)
 248
 249 static int build_cell_shell(struct cmesh *mesh, const struct cmesh *orig,
 250                 struct frac_cell *cell)
 251 {
 252         int i, j, nfaces, clipres;
 253         struct plane plane;
 254         struct poly poly, clipped, wallclipped;
 255         struct vertex *vert;
 256         int *delwall;
 257         struct plane *cplanes;
 258
 259         /* array for marking wall polygons for deletion when they get clipped entirely */
 260         delwall = alloca(cell->num_polys * sizeof *delwall);
 261         memset(delwall, 0, cell->num_polys * sizeof *delwall);
 262
 263         /* array for pre-constructing the voronoi clipping planes */
 264         cplanes = alloca(cell->num_polys * sizeof *cplanes);
 265         for(i=0; i<cell->num_polys; i++) {
 266                 poly_plane(cell->polys + i, cplanes + i);
 267                 /* flip the plane normal towards the inside of the cell, to clip everything
 268                  * outside the cell
 269                  */
 270                 cgm_vcons(&cplanes[i].norm, -cplanes[i].norm.x, -cplanes[i].norm.y, -cplanes[i].norm.z);
 271         }
 272
 273         nfaces = cmesh_poly_count(orig);
 274         for(i=0; i<nfaces; i++) {
 275                 if(mesh_poly(&poly, orig, i) == -1) {
 276                         cmesh_destroy(mesh);
 277                         return -1;
 278                 }
 279
 280                 for(j=0; j<cell->num_polys; j++) {
 281                         init_poly(&clipped);
 282                         if((clipres = clip_poly(&clipped, &poly, cplanes + j)) < 0) {
 283                                 destroy_poly(&clipped);
 284                                 break;
 285                         }
 286
 287                         /* if the plane clipped the polygon, and the two polygons intersect
 288                          * within their bounds, also clip the cell polygon by the original
 289                          * mesh polygon.
 290                          *
 291                          * TODO clipping with the polygon's plane is incorrect, and will lead
 292                          * to gaps in the cell walls when the surface is concave. We'll need
 293                          * to clip by the polygon itself, which can make the wall polygon
 294                          * concave, and will need to be split into multiple convex ones.
 295                          */
 296                         if(clipres == 0 && poly_poly(&poly, cell->polys + j)) {
 297                                 poly_plane(&poly, &plane);
 298                                 init_poly(&wallclipped);
 299
 300                                 if((clipres = clip_poly(&wallclipped, cell->polys + j, &plane)) < 0) {
 301                                         /* mark for deletion */
 302                                         delwall[j] = 1;
 303                                         destroy_poly(&wallclipped);
 304                                 } else if(clipres == 0) {
 305                                         destroy_poly(cell->polys + j);
 306                                         cell->polys[j] = wallclipped;
 307                                 } else {
 308                                         destroy_poly(&wallclipped);
 309                                 }
 310                         }
 311
 312                         destroy_poly(&poly);
 313                         poly = clipped;
 314                 }
 315
 316                 if(j < cell->num_polys) {
 317                         continue;
 318                 }
 319
 320                 vert = poly.verts + 1;
 321                 for(j=0; j<dynarr_size(poly.verts)-2; j++) {
 322                         ADD_VERTEX(mesh, poly.verts);
 323                         ADD_VERTEX(mesh, vert); vert++;
 324                         ADD_VERTEX(mesh, vert);
 325                 }
 326                 destroy_poly(&poly);
 327         }
 328
 329         /* remove deleted wall polygons */
 330         for(i=0; i<cell->num_polys; i++) {
 331                 if(delwall[i]) {
 332                         struct poly tmp = cell->polys[i];
 333                         cell->polys[i] = cell->polys[--cell->num_polys];
 334                         destroy_poly(&tmp);
 335                 }
 336         }
 337
 338         /* add wall polygons to the mesh */
 339         for(i=0; i<cell->num_polys; i++) {
 340                 vert = cell->polys[i].verts + 1;
 341                 for(j=0; j<dynarr_size(cell->polys[i].verts)-2; j++) {
 342                         ADD_VERTEX(mesh, cell->polys[i].verts);
 343                         ADD_VERTEX(mesh, vert); vert++;
 344                         ADD_VERTEX(mesh, vert);
 345                 }
 346         }
 347
 348         return 0;
 349 }
 350
 351 int frac_build_shell(struct fracture *frac)
 352 {
 353         int i, num_cells;
 354
 355         num_cells = dynarr_size(frac->cells);
 356         for(i=0; i<num_cells; i++) {
 357                 if(build_cell_shell(frac->cells[i].mesh, frac->mesh, frac->cells + i) == -1) {
 358                         return -1;
 359                 }
 360         }
 361
 362         return 0;
 363 }
 364
 365 int frac_build_walls(struct fracture *frac)
 366 {
 367         return -1;
 368 }
 369
 370 int frac_build(struct fracture *frac)
 371 {
 372         if(frac_build_cells(frac) == -1) {
 373                 return -1;
 374         }
 375         if(frac_build_shell(frac) == -1) {
 376                 return -1;
 377         }
 378         if(frac_build_walls(frac) == -1) {
 379                 return -1;
 380         }
 381         return 0;
 382 }