WebSVN – gelsvn – Diff – /branches/SPD/apps/MeshDistance/meshdist.cpp


#include <iostream>
 
#include "Util/Timer.h"
#include "Util/ArgExtracter.h"
 
#include "CGLA/Mat4x4f.h"
 
 
#include "Geometry/RGrid.h"
#include "Geometry/save_raw.h"
#include "Geometry/GridAlgorithm.h"
#include "Geometry/build_bbtree.h"
#include "Geometry/AABox.h"
 
#include "HMesh/triangulate.h"
#include "HMesh/x3d_load.h"
#include "HMesh/x3d_save.h"
 
 
 
using namespace std;
using namespace HMesh;
using namespace Geometry;
using namespace CGLA;
using namespace Util;
 
namespace
{
 
	Vec3i vol_dim(64);
	
		const Mat4x4f fit_bounding_volume(const Vec3f& p0,
																		const Vec3f& p7,
																		float buf_reg) 
	{
		Vec3f sz = p7 - p0;
		Vec3i dims = vol_dim;
		Vec3f scal_vec = (Vec3f(dims)-Vec3f(2*buf_reg+2))/sz;
		float scal = min(scal_vec[0], min(scal_vec[1],scal_vec[2]));
		
		Mat4x4f m = translation_Mat4x4f(Vec3f(0)+Vec3f(buf_reg+1));
		m *= scaling_Mat4x4f(Vec3f(scal));
		m *= translation_Mat4x4f(-p0);
		return m;
	}
 
	bool do_ray_tests = false;
	bool do_obb = true;
	bool do_aabb = false;
	bool flip_normals = false;
 
}
 
 
 
template<class BBTree>
class DistCompCache
{
	BBTree *T;
	float old_d;
	Vec3i old_p;
public:
 
	DistCompCache(BBTree* _T): T(_T), old_p(-99999) {}
 
	void operator()(const CGLA::Vec3i& pi, float& vox_val)
	{
		Vec3f p(pi);
 		if(sqr_length(pi-old_p)==1)
 			{
 				vox_val = T->compute_signed_distance(p,CGLA::sqr(1.001+fabs(old_d)));
 			}
 		else
			vox_val = T->compute_signed_distance(p);
		if(flip_normals) vox_val = -vox_val;
		old_p = pi;
		old_d = vox_val;		
	}
};
 
template<class BBTree>
class DistComp
{
	BBTree *T;
public:
 
	DistComp(BBTree* _T): T(_T) {}
 
	void operator()(const CGLA::Vec3i& pi, float& vox_val)
	{
		Vec3f p(pi);
		vox_val =  T->compute_signed_distance(p);
		if(flip_normals) vox_val = -vox_val;
	}
};
 
 
template<class BBTree>
class RayCast
{
	BBTree *T;
 
public:
	
	RayCast(BBTree* _T): T(_T) {}
 
	void operator()(const CGLA::Vec3i& pi, float& vox_val)
	{
		int n = T->intersect_cnt(Vec3f(pi), Vec3f(1,0,0));
		if(n%2==0)
			vox_val=1000;
		else
			vox_val=-1000;
	}
};
 
 
typedef RGrid<float> RGridf;
 
 
int main(int argc, char** argv)
{
	ArgExtracter ae(argc, argv);
 
	do_aabb = ae.extract("-A");
	do_ray_tests = ae.extract("-R");
	flip_normals = ae.extract("-f");
 
	Manifold m;
	x3d_load(ae.get_last_arg(), m);
 
	if(!m.is_valid())
		{
			cout << "Not a valid manifold" << endl;
			exit(0);
		}
	triangulate(m);
 
	Vec3f p0,p7;
	m.get_bbox(p0, p7);
 
	Mat4x4f T = fit_bounding_volume(p0,p7,3);
 
	for(VertexIter v = m.vertices_begin(); v != m.vertices_end(); ++v)
			v->pos = T.mul_3D_point(v->pos);
	
	
 	RGridf grid(vol_dim,FLT_MAX);
	Util::Timer tim;
 
 
	float T_build_obb=0, T_build_aabb=0, T_dist_obb=0, 
		T_dist_aabb=0, T_ray_obb=0, T_ray_aabb=0;
 
	if(do_obb)
	{
		tim.start();
		OBBTree obb_tree;
		build_OBBTree(m, obb_tree);
		T_build_obb = tim.get_secs();
		
		tim.start();
		DistCompCache<OBBTree> dist(&obb_tree);
		for_each_voxel(grid, dist);
		T_dist_obb = tim.get_secs();
		
		save_raw_float("obb_dist.raw", grid);
 
		if(do_ray_tests)
			{
				tim.start();
				RayCast<OBBTree> ray(&obb_tree);
				for_each_voxel(grid, ray);
				T_ray_obb = tim.get_secs();
				
				save_raw_float("obb_ray.raw", grid);
			}
	}
 
	if(do_aabb)
	{
		tim.start();
		AABBTree aabb_tree;
		build_AABBTree(m, aabb_tree);
		T_build_aabb = tim.get_secs();
		
		tim.start();
		DistCompCache<AABBTree> dist(&aabb_tree);
		for_each_voxel(grid, dist);
		T_dist_aabb = tim.get_secs();
		
		save_raw_float("aabb_dist.raw", grid);
 
		if(do_ray_tests)
			{
				tim.start();
				RayCast<AABBTree> ray(&aabb_tree);
				for_each_voxel(grid, ray);
				T_ray_aabb = tim.get_secs();
		
				save_raw_float("aabb_ray.raw", grid);
			}
	}
	cout.width(10);
	cout << "Poly";
	cout.width(11);
	cout <<"build_obb";
	cout.width(12);
	cout << "build_aabb";
	cout.width(10);
	cout << "dist_obb" ;
	cout.width(10);
	cout << "dist_aabb";
	cout.width(10);
	cout << "ray_obb" ;
	cout.width(10);
	cout << "ray_aabb";
	cout << endl;
 
	cout.precision(4);
	cout.width(10);
	cout << m.no_faces() << " ";
	cout.width(10);
	cout << T_build_obb;
	cout.width(12);
	cout << T_build_aabb;
	cout.width(10);
	cout << T_dist_obb;
	cout.width(10);
	cout << T_dist_aabb;
	cout.width(10);
	cout << T_ray_obb;
	cout.width(10);
	cout << T_ray_aabb;
	cout << endl;
}
 
 
 
 
 
//  	cout << "Computing distances" << endl;
// 	for(int i=0;i<triangle_vec_global.size(); ++i)
// 		{
// 			k=0;
// 			TRI = triangle_vec_global[i];
// 			for_each_voxel(grid, dist_brute);
// 		}
 
 

Rev 291	Rev 292
1	`#include <iostream>`	1	`#include <iostream>`
2		2
3	`#include "Util/Timer.h"`	3	`#include "Util/Timer.h"`
4	`#include "Util/ArgExtracter.h"`	4	`#include "Util/ArgExtracter.h"`
5		5
6	`#include "CGLA/Mat4x4f.h"`	6	`#include "CGLA/Mat4x4f.h"`
7		7
8		8
9	`#include "Geometry/RGrid.h"`	9	`#include "Geometry/RGrid.h"`
10	`#include "Geometry/save_raw.h"`	10	`#include "Geometry/save_raw.h"`
11	`#include "Geometry/GridAlgorithm.h"`	11	`#include "Geometry/GridAlgorithm.h"`
-		12	`#include "Geometry/build_bbtree.h"`
-		13	`#include "Geometry/AABox.h"`
12		14
13	`#include "HMesh/triangulate.h"`	15	`#include "HMesh/triangulate.h"`
14	`#include "HMesh/x3d_load.h"`	16	`#include "HMesh/x3d_load.h"`
15	`#include "HMesh/x3d_save.h"`	17	`#include "HMesh/x3d_save.h"`
16	`#include "build_bbtree.h"`	-
17	`#include "AABox.h"`	-
18		18
19	`using namespace std;`	19	`using namespace std;`
20	`using namespace HMesh;`	20	`using namespace HMesh;`
21	`using namespace Geometry;`	21	`using namespace Geometry;`
22	`using namespace CGLA;`	22	`using namespace CGLA;`
23	`using namespace Util;`	23	`using namespace Util;`
24		24
25	`namespace`	25	`namespace`
26	`{`	26	`{`
27		27
28	`Vec3i vol_dim(64);`	28	`Vec3i vol_dim(64);`
29		29
30	`const Mat4x4f fit_bounding_volume(const Vec3f& p0,`	30	`const Mat4x4f fit_bounding_volume(const Vec3f& p0,`
31	`const Vec3f& p7,`	31	`const Vec3f& p7,`
32	`float buf_reg)`	32	`float buf_reg)`
33	`{`	33	`{`
34	`Vec3f sz = p7 - p0;`	34	`Vec3f sz = p7 - p0;`
35	`Vec3i dims = vol_dim;`	35	`Vec3i dims = vol_dim;`
36	`Vec3f scal_vec = (Vec3f(dims)-Vec3f(2*buf_reg+2))/sz;`	36	`Vec3f scal_vec = (Vec3f(dims)-Vec3f(2*buf_reg+2))/sz;`
37	`float scal = min(scal_vec[0], min(scal_vec[1],scal_vec[2]));`	37	`float scal = min(scal_vec[0], min(scal_vec[1],scal_vec[2]));`
38		38
39	`Mat4x4f m = translation_Mat4x4f(Vec3f(0)+Vec3f(buf_reg+1));`	39	`Mat4x4f m = translation_Mat4x4f(Vec3f(0)+Vec3f(buf_reg+1));`
40	`m *= scaling_Mat4x4f(Vec3f(scal));`	40	`m *= scaling_Mat4x4f(Vec3f(scal));`
41	`m *= translation_Mat4x4f(-p0);`	41	`m *= translation_Mat4x4f(-p0);`
42	`return m;`	42	`return m;`
43	`}`	43	`}`
44		44
45	`bool do_ray_tests = false;`	45	`bool do_ray_tests = false;`
46	`bool do_obb = true;`	46	`bool do_obb = true;`
47	`bool do_aabb = false;`	47	`bool do_aabb = false;`
48	`bool flip_normals = false;`	48	`bool flip_normals = false;`
49		49
50	`}`	50	`}`
51		51
52		52
53		53
54	`template<class BBTree>`	54	`template<class BBTree>`
55	`class DistCompCache`	55	`class DistCompCache`
56	`{`	56	`{`
57	`BBTree *T;`	57	`BBTree *T;`
58	`float old_d;`	58	`float old_d;`
59	`Vec3i old_p;`	59	`Vec3i old_p;`
60	`public:`	60	`public:`
61		61
62	`DistCompCache(BBTree* _T): T(_T), old_p(-99999) {}`	62	`DistCompCache(BBTree* _T): T(_T), old_p(-99999) {}`
63		63
64	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`	64	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`
65	`{`	65	`{`
66	`Vec3f p(pi);`	66	`Vec3f p(pi);`
67	`if(sqr_length(pi-old_p)==1)`	67	`if(sqr_length(pi-old_p)==1)`
68	`{`	68	`{`
69	`vox_val = T->compute_signed_distance(p,CGLA::sqr(1.001+fabs(old_d)));`	69	`vox_val = T->compute_signed_distance(p,CGLA::sqr(1.001+fabs(old_d)));`
70	`}`	70	`}`
71	`else`	71	`else`
72	`vox_val = T->compute_signed_distance(p);`	72	`vox_val = T->compute_signed_distance(p);`
73	`if(flip_normals) vox_val = -vox_val;`	73	`if(flip_normals) vox_val = -vox_val;`
74	`old_p = pi;`	74	`old_p = pi;`
75	`old_d = vox_val;`	75	`old_d = vox_val;`
76	`}`	76	`}`
77	`};`	77	`};`
78		78
79	`template<class BBTree>`	79	`template<class BBTree>`
80	`class DistComp`	80	`class DistComp`
81	`{`	81	`{`
82	`BBTree *T;`	82	`BBTree *T;`
83	`public:`	83	`public:`
84		84
85	`DistComp(BBTree* _T): T(_T) {}`	85	`DistComp(BBTree* _T): T(_T) {}`
86		86
87	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`	87	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`
88	`{`	88	`{`
89	`Vec3f p(pi);`	89	`Vec3f p(pi);`
90	`vox_val = T->compute_signed_distance(p);`	90	`vox_val = T->compute_signed_distance(p);`
91	`if(flip_normals) vox_val = -vox_val;`	91	`if(flip_normals) vox_val = -vox_val;`
92	`}`	92	`}`
93	`};`	93	`};`
94		94
95		95
96	`template<class BBTree>`	96	`template<class BBTree>`
97	`class RayCast`	97	`class RayCast`
98	`{`	98	`{`
99	`BBTree *T;`	99	`BBTree *T;`
100		100
101	`public:`	101	`public:`
102		102
103	`RayCast(BBTree* _T): T(_T) {}`	103	`RayCast(BBTree* _T): T(_T) {}`
104		104
105	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`	105	`void operator()(const CGLA::Vec3i& pi, float& vox_val)`
106	`{`	106	`{`
107	`int n = T->intersect_cnt(Vec3f(pi), Vec3f(1,0,0));`	107	`int n = T->intersect_cnt(Vec3f(pi), Vec3f(1,0,0));`
108	`if(n%2==0)`	108	`if(n%2==0)`
109	`vox_val=1000;`	109	`vox_val=1000;`
110	`else`	110	`else`
111	`vox_val=-1000;`	111	`vox_val=-1000;`
112	`}`	112	`}`
113	`};`	113	`};`
114		114
115		115
116	`typedef RGrid<float> RGridf;`	116	`typedef RGrid<float> RGridf;`
117		117
118		118
119	`int main(int argc, char** argv)`	119	`int main(int argc, char** argv)`
120	`{`	120	`{`
121	`ArgExtracter ae(argc, argv);`	121	`ArgExtracter ae(argc, argv);`
122		122
123	`do_aabb = ae.extract("-A");`	123	`do_aabb = ae.extract("-A");`
124	`do_ray_tests = ae.extract("-R");`	124	`do_ray_tests = ae.extract("-R");`
125	`flip_normals = ae.extract("-f");`	125	`flip_normals = ae.extract("-f");`
126		126
127	`Manifold m;`	127	`Manifold m;`
128	`x3d_load(ae.get_last_arg(), m);`	128	`x3d_load(ae.get_last_arg(), m);`
129		129
130	`if(!m.is_valid())`	130	`if(!m.is_valid())`
131	`{`	131	`{`
132	`cout << "Not a valid manifold" << endl;`	132	`cout << "Not a valid manifold" << endl;`
133	`exit(0);`	133	`exit(0);`
134	`}`	134	`}`
135	`triangulate(m);`	135	`triangulate(m);`
136		136
137	`Vec3f p0,p7;`	137	`Vec3f p0,p7;`
138	`m.get_bbox(p0, p7);`	138	`m.get_bbox(p0, p7);`
139		139
140	`Mat4x4f T = fit_bounding_volume(p0,p7,3);`	140	`Mat4x4f T = fit_bounding_volume(p0,p7,3);`
141		141
142	`for(VertexIter v = m.vertices_begin(); v != m.vertices_end(); ++v)`	142	`for(VertexIter v = m.vertices_begin(); v != m.vertices_end(); ++v)`
143	`v->pos = T.mul_3D_point(v->pos);`	143	`v->pos = T.mul_3D_point(v->pos);`
144		144
145		145
146	`RGridf grid(vol_dim,FLT_MAX);`	146	`RGridf grid(vol_dim,FLT_MAX);`
147	`Util::Timer tim;`	147	`Util::Timer tim;`
148		148
149		149
150	`float T_build_obb=0, T_build_aabb=0, T_dist_obb=0,`	150	`float T_build_obb=0, T_build_aabb=0, T_dist_obb=0,`
151	`T_dist_aabb=0, T_ray_obb=0, T_ray_aabb=0;`	151	`T_dist_aabb=0, T_ray_obb=0, T_ray_aabb=0;`
152		152
153	`if(do_obb)`	153	`if(do_obb)`
154	`{`	154	`{`
155	`tim.start();`	155	`tim.start();`
156	`OBBTree obb_tree;`	156	`OBBTree obb_tree;`
157	`build_OBBTree(m, obb_tree);`	157	`build_OBBTree(m, obb_tree);`
158	`T_build_obb = tim.get_secs();`	158	`T_build_obb = tim.get_secs();`
159		159
160	`tim.start();`	160	`tim.start();`
161	`DistCompCache<OBBTree> dist(&obb_tree);`	161	`DistCompCache<OBBTree> dist(&obb_tree);`
162	`for_each_voxel(grid, dist);`	162	`for_each_voxel(grid, dist);`
163	`T_dist_obb = tim.get_secs();`	163	`T_dist_obb = tim.get_secs();`
164		164
165	`save_raw_float("obb_dist.raw", grid);`	165	`save_raw_float("obb_dist.raw", grid);`
166		166
167	`if(do_ray_tests)`	167	`if(do_ray_tests)`
168	`{`	168	`{`
169	`tim.start();`	169	`tim.start();`
170	`RayCast<OBBTree> ray(&obb_tree);`	170	`RayCast<OBBTree> ray(&obb_tree);`
171	`for_each_voxel(grid, ray);`	171	`for_each_voxel(grid, ray);`
172	`T_ray_obb = tim.get_secs();`	172	`T_ray_obb = tim.get_secs();`
173		173
174	`save_raw_float("obb_ray.raw", grid);`	174	`save_raw_float("obb_ray.raw", grid);`
175	`}`	175	`}`
176	`}`	176	`}`
177		177
178	`if(do_aabb)`	178	`if(do_aabb)`
179	`{`	179	`{`
180	`tim.start();`	180	`tim.start();`
181	`AABBTree aabb_tree;`	181	`AABBTree aabb_tree;`
182	`build_AABBTree(m, aabb_tree);`	182	`build_AABBTree(m, aabb_tree);`
183	`T_build_aabb = tim.get_secs();`	183	`T_build_aabb = tim.get_secs();`
184		184
185	`tim.start();`	185	`tim.start();`
186	`DistCompCache<AABBTree> dist(&aabb_tree);`	186	`DistCompCache<AABBTree> dist(&aabb_tree);`
187	`for_each_voxel(grid, dist);`	187	`for_each_voxel(grid, dist);`
188	`T_dist_aabb = tim.get_secs();`	188	`T_dist_aabb = tim.get_secs();`
189		189
190	`save_raw_float("aabb_dist.raw", grid);`	190	`save_raw_float("aabb_dist.raw", grid);`
191		191
192	`if(do_ray_tests)`	192	`if(do_ray_tests)`
193	`{`	193	`{`
194	`tim.start();`	194	`tim.start();`
195	`RayCast<AABBTree> ray(&aabb_tree);`	195	`RayCast<AABBTree> ray(&aabb_tree);`
196	`for_each_voxel(grid, ray);`	196	`for_each_voxel(grid, ray);`
197	`T_ray_aabb = tim.get_secs();`	197	`T_ray_aabb = tim.get_secs();`
198		198
199	`save_raw_float("aabb_ray.raw", grid);`	199	`save_raw_float("aabb_ray.raw", grid);`
200	`}`	200	`}`
201	`}`	201	`}`
202	`cout.width(10);`	202	`cout.width(10);`
203	`cout << "Poly";`	203	`cout << "Poly";`
204	`cout.width(11);`	204	`cout.width(11);`
205	`cout <<"build_obb";`	205	`cout <<"build_obb";`
206	`cout.width(12);`	206	`cout.width(12);`
207	`cout << "build_aabb";`	207	`cout << "build_aabb";`
208	`cout.width(10);`	208	`cout.width(10);`
209	`cout << "dist_obb" ;`	209	`cout << "dist_obb" ;`
210	`cout.width(10);`	210	`cout.width(10);`
211	`cout << "dist_aabb";`	211	`cout << "dist_aabb";`
212	`cout.width(10);`	212	`cout.width(10);`
213	`cout << "ray_obb" ;`	213	`cout << "ray_obb" ;`
214	`cout.width(10);`	214	`cout.width(10);`
215	`cout << "ray_aabb";`	215	`cout << "ray_aabb";`
216	`cout << endl;`	216	`cout << endl;`
217		217
218	`cout.precision(4);`	218	`cout.precision(4);`
219	`cout.width(10);`	219	`cout.width(10);`
220	`cout << m.no_faces() << " ";`	220	`cout << m.no_faces() << " ";`
221	`cout.width(10);`	221	`cout.width(10);`
222	`cout << T_build_obb;`	222	`cout << T_build_obb;`
223	`cout.width(12);`	223	`cout.width(12);`
224	`cout << T_build_aabb;`	224	`cout << T_build_aabb;`
225	`cout.width(10);`	225	`cout.width(10);`
226	`cout << T_dist_obb;`	226	`cout << T_dist_obb;`
227	`cout.width(10);`	227	`cout.width(10);`
228	`cout << T_dist_aabb;`	228	`cout << T_dist_aabb;`
229	`cout.width(10);`	229	`cout.width(10);`
230	`cout << T_ray_obb;`	230	`cout << T_ray_obb;`
231	`cout.width(10);`	231	`cout.width(10);`
232	`cout << T_ray_aabb;`	232	`cout << T_ray_aabb;`
233	`cout << endl;`	233	`cout << endl;`
234	`}`	234	`}`
235		235
236		236
237		237
238		238
239		239
240	`// cout << "Computing distances" << endl;`	240	`// cout << "Computing distances" << endl;`
241	`// for(int i=0;i<triangle_vec_global.size(); ++i)`	241	`// for(int i=0;i<triangle_vec_global.size(); ++i)`
242	`// {`	242	`// {`
243	`// k=0;`	243	`// k=0;`
244	`// TRI = triangle_vec_global[i];`	244	`// TRI = triangle_vec_global[i];`
245	`// for_each_voxel(grid, dist_brute);`	245	`// for_each_voxel(grid, dist_brute);`
246	`// }`	246	`// }`
247		247
248		248

Subversion Repositories gelsvn

(root)/branches/SPD/apps/MeshDistance/meshdist.cpp – Rev 291 → 292