1 /*************************************************************************
2 * *
3 * Vega FEM Simulation Library Version 3.1 *
4 * *
5 * "distance field" library , Copyright (C) 2007 CMU, 2016 USC *
6 * All rights reserved. *
7 * *
8 * Code authors: Jernej Barbic, Hongyi Xu, Yijing Li *
9 * http://www.jernejbarbic.com/code *
10 * *
11 * Research: Jernej Barbic, Hongyi Xu, Doug L. James *
12 * *
13 * Funding: National Science Foundation, Link Foundation, *
14 * Zumberge Research and Innovation Fund at USC *
15 * *
16 * This library is free software; you can redistribute it and/or *
17 * modify it under the terms of the BSD-style license that is *
18 * included with this library in the file LICENSE.txt *
19 * *
20 * This library is distributed in the hope that it will be useful, *
21 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
22 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the file *
23 * LICENSE.TXT for more details. *
24 * *
25 *************************************************************************/
26
27 /*
28 This code generates a distance filed, either signed or unsigned,
29 to the given triangle mesh.
30
31 The distance field can be loaded/saved to a file.
32 You can also lookup the field (once computed or loaded) at arbitrary
33 locations inside the field, using trilinear interpolation.
34
35 Input mesh need not have triangle faces (e.g., can have quads;
36 they will be triangulated).
37
38 For signed field generation, the input mesh must be a closed manifold mesh.
39
40 Input mesh must be given in the .obj format:
41 http://www.royriggs.com/obj.html
42
43 By default, the bounding box will be a cube, obtained by fitting the
44 smallest cube to the geometry, and then expanded (scaled) from its
45 center by a factor of 1.5. You can provide your own bounding boxes.
46 However, note: (1) the provided bounding box must cover all the geometry,
47 and (2) bounding boxes that are not cubes were not (yet) tested.
48
49 The bounding box will be divided into a "resolution" number of
50 cubes ("voxels") along each axis. The distance field will be computed at the
51 vertices of these voxels. So, if resolution is 256, the bounding box
52 will get divided into 256 x 256 x 256 voxels, and the distance field
53 will be computed on the resulting 257 x 257 x 257 grid of voxel vertices.
54 Note that when indexing voxels, the indices (i,j,k) will run from 0 to 255
55 inclusive, whereas when indexing voxel vertices (also called "grid vertices"),
56 they will run from 0 to 256 inclusive.
57
58 Distance field data is stored at voxel vertices.
59 In memory, distance field value at voxel vertex (i,j,k) is stored
60 at location k * (resolutionX+1)*(resolutionY+1) + j * (resolutionX+1) + i .
61
62 Internally, the code builds an octree on top of the triangle mesh.
63 There are two parameters that control this process (you can keep them
64 at default values, which worked well in practice for us) :
65 the max depth of the octree is "maxDepth", and
66 the max number of triangles intersecting an octree cell is "maxTriCount".
67 Note: once max depth level is reached, the maxTriCount bound is not imposed any more.
68 */
69
70 #include <iostream>
71 #include "objMeshOrientable.h"
72 #include "objMeshOctree.h"
73 #include "distanceFieldBase.h"
74 using namespace std;
75
DistanceFieldBase()76 DistanceFieldBase::DistanceFieldBase()
77 {
78 resolutionX = 1;
79 resolutionY = 1;
80 resolutionZ = 1;
81 setAutomaticBoundingBox();
82 bmin_ = Vec3d(0,0,0);
83 bmax_ = Vec3d(1,1,1);
84 setGridParameters();
85 bboxComputed = false;
86 }
87
DistanceFieldBase(int resolutionX_,int resolutionY_,int resolutionZ_)88 DistanceFieldBase::DistanceFieldBase(int resolutionX_, int resolutionY_, int resolutionZ_) : resolutionX(resolutionX_), resolutionY(resolutionY_), resolutionZ(resolutionZ_)
89 {
90 setAutomaticBoundingBox();
91 bmin_ = Vec3d(0,0,0);
92 bmax_ = Vec3d(1,1,1);
93 setGridParameters();
94 bboxComputed = false;
95 }
96
computeBoundingBox(ObjMesh * objMesh,int resolutionX,int resolutionY,int resolutionZ)97 void DistanceFieldBase::computeBoundingBox(ObjMesh * objMesh, int resolutionX, int resolutionY, int resolutionZ)
98 {
99 Vec3d bminTemp, bmaxTemp;
100 objMesh->getBoundingBox(1.0, &bminTemp, &bmaxTemp);
101
102 // set aspect ratio that corresponds to the resolutions
103 Vec3d bcenterTemp = 0.5 * (bminTemp + bmaxTemp);
104
105 cout << "Tight bounding box:" << endl << " " << bminTemp << endl << " " << bmaxTemp << endl;
106
107 Vec3d sideTemp = bmaxTemp - bminTemp;
108 if (sideTemp[0] / resolutionX < sideTemp[1] / resolutionY)
109 {
110 // increase x
111 sideTemp[0] = sideTemp[1] / resolutionY * resolutionX;
112 }
113 else
114 {
115 // increase y
116 sideTemp[1] = sideTemp[0] / resolutionX * resolutionY;
117 }
118
119 // now x,y are ok, must adjust z
120 if (sideTemp[1] / resolutionY < sideTemp[2] / resolutionZ)
121 {
122 // increase x and y
123 double factor = (sideTemp[2] / resolutionZ * resolutionY) / sideTemp[1];
124 sideTemp[1] *= factor;
125 sideTemp[0] *= factor;
126 }
127 else
128 {
129 // increase z
130 sideTemp[2] = sideTemp[1] / resolutionY * resolutionZ;
131 }
132
133 BoundingBox bbox(bminTemp, bmaxTemp);
134 bbox.setbmin(bcenterTemp - 0.5 * sideTemp);
135 bbox.setbmax(bcenterTemp + 0.5 * sideTemp);
136 if (allBoxSidesEqual)
137 bbox.regularize();
138 bbox.expand(expansionRatio);
139
140 bmin_ = bbox.bmin();
141 bmax_ = bbox.bmax();
142
143 bboxComputed = true;
144 }
145
146