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: Hongyi Xu, Jernej Barbic *
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 Narrowband distance field computation.
29 */
30
31 #ifndef _DISTANCEFIELDNARROWBAND_H_
32 #define _DISTANCEFIELDNARROWBAND_H_
33
34 #include "objMesh.h"
35 #include "triple.h"
36 #include "distanceFieldBase.h"
37 #include "vegalong.h"
38 #include <map>
39 #include <float.h>
40
41 class DistanceFieldNarrowBand : public DistanceFieldBase
42 {
43 public:
44
45 DistanceFieldNarrowBand();
46 virtual ~DistanceFieldNarrowBand();
47
48 // computes unsigned distance field in a narrow band
49 // filename is the obj filename
50 // bandWidth is given in absolute units
51 virtual int computeUnsignedField(ObjMesh * objMesh, int resolutionX, int resolutionY, int resolutionZ, double bandWidth, int maxTriCount=15, int maxDepth=10);
52
53 // computes signed distance field in a narrow band
54 // filename is the obj filename
55 // offset is given in absolute units
56 virtual int computeSignedField(ObjMesh * objMesh, int resolutionX, int resolutionY, int resolutionZ, double bandWidth, int maxTriCount=15, int maxDepth=10);
57
58 virtual void offsetDistanceField(double sigma); // add sigma to all the distance field values
59 // computes the signed distance field, in the interior region (bounded by the zero isosurface), in a narrow band
60 // must call computeUnsignedField using the same objMesh and same resolution first
61 // also, you should first offset the unsigned distance field by a proper 'sigma' (typically, by passing -sigma, where sigma > 0)
62 virtual int computeInteriorSignedField(ObjMesh * objMesh, int resolutionX, int resolutionY, int resolutionZ, double bandWidth, int maxTriCount=15, int maxDepth=10);
63
64 // loads a previously computed distance field from a disk file
65 virtual int load(const std::string& filename); // returns 0 on success
66
67 // saves the current distance field to a disk file (e.g. after computing it once, for later fast reloading)
68 virtual int save(const std::string& filename, bool doublePrecision); // saves in a compressed format (saves the computed grid points only)
69 virtual int saveToText(const std::string& filename);
70 int saveToDistanceField(const std::string& filename, bool doublePrecision); // saves in full format (FLT_MAX or -FLT_MAX are stored for uncomputed grid points)
71
72 // return distance field value at grid vertex (i,j,k)
73 // each of i,j,k must be an integer from {0, ..., resolution{X,Y,Z}}
74 virtual inline float distance(int i, int j, int k);
75 // computes distance and gradient at arbitrary position
76 virtual float distance(Vec3d pos, int constrainToBox=0);
77 // alters the distance at a particular grid vertex (i,j,k)
78 virtual inline void setDistance(int i, int j, int k, float value);
79
80 virtual Vec3d gradient(const Vec3d & pos);
81
82 virtual bool sanityCheck(); // checks if distance for any two adjacent voxels is less than voxel grid spacing apart (which it must be by triangle inequality, for both signed and unsigned fields)
83
84 virtual float maxValue();
85 virtual float minValue();
86 virtual void maxMinValue(float* maxValue, float* minValue);
87
88 virtual float maxAbsValue();
89 virtual float maxAbsValue(float threshold); // only abs values up to threshold
90 virtual float maxNonInftyAbsValue();
91
92 virtual int breadthFirstTraversalSigned(void * objMeshOctree, float offset, int zLo, int zHi, int asterisk=0);
93 virtual int breadthFirstTraversalUnsigned(void * objMeshOctree, float offset, int zLo, int zHi, int asterisk=0);
94 virtual int breadthFirstTraversalInteriorSigned(void * objMeshOctree, float offset, int zLo, int zHi, int asterisk=0);
95
96 void freeMemory();
97
98 void findSurfaceGridPoints(ObjMesh* objMeshIn);
99
100 typedef triple<int,int,int> gridPoint;
getDistanceData()101 std::map<gridPoint, float> * getDistanceData() {return &distanceData;};
102
103 protected:
104 int maxTriCount;
105 int maxDepth;
106
107 void finalizeGridPointStatus();
108 enum { COMPUTED, EXTERIOR_UNCOMPUTED, INTERIOR_UNCOMPUTED };
109 char * gridPointStatus; // needed for the sign on uncomputed grid points; unsigned field only has two values: COMPUTED, EXTERIOR_UNCOMPUTED
110
111 int signFieldFlooded;
112 //float * fieldData;
113 std::map<gridPoint, float> distanceData;
114
115 vegalong GetFilesize(const char *filename);
116
117 std::vector<vegalong> surfaceGridPoints;
118 };
119
distance(int i,int j,int k)120 inline float DistanceFieldNarrowBand::distance(int i, int j, int k)
121 {
122 vegalong index = (k * (resolutionY+1) + j) * (resolutionX + 1) + i;
123
124 std::map<gridPoint, float>::iterator it = distanceData.find(gridPoint(i,j,k));
125
126 if (it == distanceData.end())
127 {
128 if (gridPointStatus[index] == DistanceFieldNarrowBand::INTERIOR_UNCOMPUTED)
129 return -FLT_MAX;
130 else
131 return FLT_MAX;
132 }
133 else
134 return it->second;
135 }
136
setDistance(int i,int j,int k,float value)137 inline void DistanceFieldNarrowBand::setDistance(int i, int j, int k, float value)
138 {
139 std::map<gridPoint, float>::iterator it = distanceData.find(gridPoint(i,j,k));
140
141 if (it == distanceData.end())
142 distanceData.insert(std::pair<gridPoint, float>(gridPoint(i,j,k), value));
143 else
144 it->second = value;
145 }
146
147 #endif
148
149