Subversion Repositories gelsvn

Rev

Rev 630 | Go to most recent revision | Only display areas with differences | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

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