WebSVN – gelsvn – Diff – /trunk/src/Geometry/tessellate.cpp


/* ----------------------------------------------------------------------- *
 * This file is part of GEL, http://www.imm.dtu.dk/GEL
 * Copyright (C) the authors and DTU Informatics
 * For license and list of authors, see ../../doc/intro.pdf
 * ----------------------------------------------------------------------- */
 
#include <map>
#include <algorithm>
#include <queue>
 
#include "Util/Grid2D.h"
#include "CGLA/Vec3f.h"
#include "CGLA/Vec2f.h"
 
#include "tessellate.h"
 
using namespace std;
using namespace CGLA;
using namespace Util;
 
namespace Geometry
{
	
	float MAX_ERR = 5;
	float MAX_DIST = 5;
	int  ADAPTIVE = false;
	
#define mymax(x,y) ((x)>(y)?(x):(y))
	
	// Test surfaces ---
	
	
	// ----------- CONSTANTS AND GLOBALS ------------
	
	namespace
	{
		float maximum_edge_error=0;
		
		IndexedFaceSet* face_set_ptr = 0;
		ParSurf* the_surf = 0;
		
		class TreeNode;
		typedef TreeNode* TreeNodePtr;
		
		vector<TreeNode*> root_nodes;
		
		class Point 
		{
			friend class Edge;
			friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);
			friend float compute_mid_pt_err(const Point&, const Point&, Point&);
			
			
			Vec2f ppos;
			Vec3f pos;
			mutable int vertex_idx;
			
public:
				
				Point(const Vec2f& _ppos, const Vec3f& _pos): 
				ppos(_ppos), pos(_pos), vertex_idx(-1)
			{
			}
			
			Point() {}
			
			void create_vertex() const
			{
				if(vertex_idx == -1) 
					vertex_idx = face_set_ptr->add_vertex(pos);
			}
		};
		
		
		
		inline const Point create_point(float u, float v)
		{
			Vec2f ppos(u,v);
			Vec3f pos = (*the_surf)(u,v);
			return Point(ppos, pos);
		}
		
		float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node);
		
		class TreeNode
		{
			friend class Edge;
			friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);
			
			const Point p;
			const float error;
			
			const TreeNodePtr left;
			const TreeNodePtr right;
			
public:
				
				TreeNode(const Point& _p, float _error,
						 TreeNodePtr _left, TreeNodePtr _right):
				p(_p), error(_error), left(_left), right(_right)  {}
			
			~TreeNode() 
			{
				if(left)
					delete left; 
				
				if(right)
					delete right;
			}
			
		};
		
		class Edge
		{
			friend Edge create_edge(const Point&, const Point&, bool);
			
			Edge(const Point& _left_point, 
				 const Point& _right_point, 
				 TreeNodePtr _root,
				 bool _cw): 
				left_point(_left_point), right_point(_right_point), root(_root), cw(_cw)
			{
			}
			
			
			Point left_point;
			Point right_point;
			TreeNodePtr root;
			bool cw;
			
public:
				
				Edge() {}
			
			float get_error() const
			{
				if(root)
					return root->error;
				return 0.0f;
			}
			
			float get_length() const
			{
				return (left_point.pos-right_point.pos).length();
			}
			
			bool is_simple() const 
			{
				if(root->left==0)
				{
					assert(root->left==0);
					assert(root->right==0);
					return true;
				}
				return false;
			}
			
			const Edge get_sub_edge(int edge_no) const
			{
				assert(edge_no == 0 || edge_no == 1);
				assert(!is_simple());
				
				if(!cw) edge_no = 1-edge_no;
				
				if(edge_no ==0)
					return Edge(left_point, root->p, root->left, cw);
				else
					return Edge(root->p, right_point, root->right, cw);
			}
			
			const Edge get_opp_edge() const
			{
				return Edge(left_point, right_point, root, !cw);
			}
			
			
			const Point& get_split_point() const
			{
				assert(!is_simple());
				return root->p;
			}
			
			const Point& get_point(int i) const
			{
				assert(i==0 || i==1);
				if(!cw) i=1-i;
				if(i == 0) 
					return left_point;
				return right_point;
			}
			
			int get_point_idx(int i) const
			{
				assert(i==0 || i==1);
				if(!cw) i=1-i;
				if(i == 0) 
					return left_point.vertex_idx;
				return right_point.vertex_idx;
			}
			
		};
		
		
		float compute_mid_pt_err(const Point&, const Point&, Point&);
		
		Edge create_edge(const Point& left_point, const Point& right_point, 
						 bool cw=true)
		{
			TreeNodePtr root=0;
			if(ADAPTIVE)
				subdiv(0, left_point,right_point,root);
			else
			{
				Point mid_point;
				float err = compute_mid_pt_err(left_point, right_point, mid_point);
				root = new TreeNode(mid_point,err,0,0);
			}
			if(root)
				root_nodes.push_back(root);
			
			return Edge(left_point, right_point, root, cw);
		}
		
		class Triangle
		{
			Edge a,b,c;
public:
			Triangle(const Edge& _a, const Edge& _b, const Edge& _c):
			a(_a),  b(_b), c(_c) {}
			
			const Edge& operator[](int i) const
		{
			switch(i)
			{
				case 0: return a;
				case 1: return b;
				case 2: return c;
			}
			return c;
		}
			
		};
		
		float compute_mid_pt_err(const Point& lp, const Point& rp, 
								 Point& mid_point)
		{
			Vec3f m = (lp.pos - rp.pos)/2.0f + rp.pos;
			Vec2f uvm = (lp.ppos -rp.ppos)/2.0f + rp.ppos;
			
			mid_point = create_point(uvm[0], uvm[1]);
			Vec3f diff = mid_point.pos - m;
			return diff.length();
		}
		
		float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node)
		{
			Point mid_point;
			float err = compute_mid_pt_err(lp,rp,mid_point);
			
//			Vec3f dist = lp.pos-rp.pos;
			if (level>10) 
			{
				node = new TreeNode(mid_point,err,0,0);
				return err;
			}
			
			TreeNodePtr lnode;
			float lerr = subdiv(level+1,lp, mid_point, lnode);
			
			TreeNodePtr rnode;
			float rerr = subdiv(level+1,mid_point, rp, rnode);
			
			float max_err = mymax(err, mymax(lerr,rerr));
			
			if(max_err < MAX_ERR)
			{
				if(lnode) delete lnode;
				if(rnode) delete rnode;
				
				node = new TreeNode(mid_point, max_err, 0,0);
			}
			else
			{
				mid_point.create_vertex();
				node = new TreeNode(mid_point, max_err, lnode,rnode);
			}
			
			return max_err;
		}
		
		
		
		void split(const Triangle orig_tri)
		{
			queue<Triangle> tri_queue;
			tri_queue.push(orig_tri);
			
			
			while(!tri_queue.empty())
			{
				const Triangle tri = tri_queue.front();
				tri_queue.pop();
				
				int no_split_edges = 0;
				
				for(int i=0;i<3;++i)
					if(!tri[i].is_simple())
						++no_split_edges;
				
				if(no_split_edges==0 /* || level > 2000*/) 
				{
					int a = tri[0].get_point_idx(0);
					int b = tri[1].get_point_idx(0);
					int c = tri[2].get_point_idx(0);
					face_set_ptr->add_face(Vec3i(a,b,c));
					maximum_edge_error = mymax(mymax(tri[0].get_error(),
													 tri[1].get_error()),
											   tri[2].get_error());
				}
				else if(no_split_edges==3)
				{
					const Point& p0 = tri[0].get_point(0);
					const Point& p1 = tri[1].get_point(0);
					const Point& p2 = tri[2].get_point(0);
					
					const Point sp[3] = {tri[0].get_split_point(),
						tri[1].get_split_point(),
						tri[2].get_split_point()};
					
					Edge sp_edg[3] = {create_edge(sp[0],p2),
						create_edge(sp[1],p0),
						create_edge(sp[2],p1)};
					
					int i_min=0;
					float l_min = sp_edg[0].get_length();
					int i;
					for(i=1;i<3;++i)
					{
						float len = sp_edg[i].get_length();
						if(len<l_min) 
						{
							l_min = len;
							i_min = i;
						}
					}
					
					i = i_min;
					int j = (i+1)%3;
					int k = (i+2)%3;
					
					Edge spi_spj = create_edge(sp[i],sp[j]);
					Edge spi_spk = create_edge(sp[i],sp[k]);
					
					Triangle tri0(tri[i].get_sub_edge(1),
								  tri[j].get_sub_edge(0),
								  spi_spj.get_opp_edge());
					Triangle tri1(spi_spj,
								  tri[j].get_sub_edge(1),
								  sp_edg[i].get_opp_edge());
					Triangle tri2(sp_edg[i],
								  tri[k].get_sub_edge(0),
								  spi_spk.get_opp_edge());
					Triangle tri3(spi_spk,
								  tri[k].get_sub_edge(1),
								  tri[i].get_sub_edge(0));
					
					tri_queue.push(tri0);
					tri_queue.push(tri1);
					tri_queue.push(tri2);
					tri_queue.push(tri3);
				}
				else if(no_split_edges==2)
				{
					int i;
					if(tri[2].is_simple()) i = 0;
					else if(tri[1].is_simple()) i=2;
					else i = 1;
					
					int j = (i+1)%3;
					int k = (i+2)%3;
					
					const Point& spi = tri[i].get_split_point();
					const Point& spj = tri[j].get_split_point();
					
					const Point& pi = tri[i].get_point(0);
					const Point& pk = tri[k].get_point(0);
					
					Edge spi_spj = create_edge(spi, spj);
					Edge pi_spj = create_edge(pi, spj);
					Edge pk_spi = create_edge(pk, spi);
					
					Triangle tri0(tri[i].get_sub_edge(1),
								  tri[j].get_sub_edge(0),
								  spi_spj.get_opp_edge());
					
					
					if(pi_spj.get_length() < pk_spi.get_length())
					{
						Triangle tri1(spi_spj,
									  pi_spj.get_opp_edge(),
									  tri[i].get_sub_edge(0));
						Triangle tri2(pi_spj,
									  tri[j].get_sub_edge(1),
									  tri[k]);
						tri_queue.push(tri1);					
						tri_queue.push(tri2);	
					}
					else
					{
						Triangle tri1(tri[i].get_sub_edge(0),
									  pk_spi.get_opp_edge(),
									  tri[k]);
						Triangle tri2(spi_spj,
									  tri[j].get_sub_edge(1),
									  pk_spi);
						tri_queue.push(tri1);					
						tri_queue.push(tri2);
					}
					tri_queue.push(tri0);
				}
				else if(no_split_edges==1)
				{
					int i;
					if(!tri[0].is_simple()) i=0;
					else if(!tri[1].is_simple()) i=1;
					else i= 2;
					
					int j = (i+1)%3;
					int k = (i+2)%3;
					
					const Point& pk = tri[k].get_point(0);
					const Point& spi = tri[i].get_split_point();
					
					Edge spi_pk = create_edge(spi, pk);
					
					Triangle tri0(tri[i].get_sub_edge(1),
								  tri[j],
								  spi_pk.get_opp_edge());
					
					Triangle tri1(tri[i].get_sub_edge(0),
								  spi_pk,
								  tri[k]);
					
					tri_queue.push(tri0);
					tri_queue.push(tri1);
					
				}
			}
		}
	}
	
	
	void tessellate(IndexedFaceSet& face_set, ParSurf& s, Grid2D<Vec3f>& inigrid)
	{	
		face_set_ptr = &face_set;
		int i;
		int n=inigrid.get_xdim()-1;
		int m=inigrid.get_ydim()-1;
		the_surf = &s;
		Grid2D<Point> points(n+1,m+1);
		for(i=0;i<=n;++i)
			for(int j=0;j<=m;++j)
			{
				Vec3f p = inigrid(i,j);
				points(i,j) = create_point(p[0], p[1]);
				points(i,j).create_vertex();
			}
				
				
				Grid2D<Edge> h_edges(n,m+1);
		Grid2D<Edge> v_edges(n+1,m);
		for(i=0;i<=n;++i)
			for(int j=0;j<=m;++j)
			{
				if(i<n)
					h_edges(i,j) = create_edge(points(i,j), points(i+1,j));
				if(j<m)
					v_edges(i,j) = create_edge(points(i,j), points(i,j+1));
			}
				
				for(i=0;i<n;++i)
					for(int j=0;j<m;++j)
					{
						Edge diag = create_edge(points(i,j), points(i+1,j+1));
						Triangle tri0(h_edges(i,j), v_edges(i+1,j), diag.get_opp_edge());
						Triangle tri1(diag, h_edges(i,j+1).get_opp_edge(), 
									  v_edges(i,j).get_opp_edge());
						split(tri0);
						split(tri1);
					}
						for(i=0;i<static_cast<int>(root_nodes.size());++i)
						{
							if(root_nodes[i])
							{
								delete root_nodes[i];
							}
						}
						root_nodes.clear();
		maximum_edge_error=0;
		the_surf = 0;
	}
	
	
	void tessellate(IndexedFaceSet& face_set, ParSurf& s,
					float u_min, float u_max, float v_min, float v_max, 
					int n, int m)
	{	
		face_set_ptr = &face_set;
		int i;
		the_surf = &s;
		Grid2D<Point> points(n+1,m+1);
		float step_x = (u_max-u_min)/float(n);
		float step_y = (v_max-v_min)/float(m);
		for(i=0;i<=n;++i)
			for(int j=0;j<=m;++j)
			{
				points(i,j) = create_point(u_min+i*step_x,v_min+j*step_y);
				points(i,j).create_vertex();
			}
				
				
				Grid2D<Edge> h_edges(n,m+1);
		Grid2D<Edge> v_edges(n+1,m);
		for(i=0;i<=n;++i)
			for(int j=0;j<=m;++j)
			{
				if(i<n)
					h_edges(i,j) = create_edge(points(i,j), points(i+1,j));
				if(j<m)
					v_edges(i,j) = create_edge(points(i,j), points(i,j+1));
			}
				
				for(i=0;i<n;++i)
					for(int j=0;j<m;++j)
					{
						Edge diag = create_edge(points(i,j), points(i+1,j+1));
						Triangle tri0(h_edges(i,j), v_edges(i+1,j), diag.get_opp_edge());
						Triangle tri1(diag, h_edges(i,j+1).get_opp_edge(), 
									  v_edges(i,j).get_opp_edge());
						split(tri0);
						split(tri1);
					}
						for(i=0;i<static_cast<int>(root_nodes.size());++i)
						{
							if(root_nodes[i])
							{
								delete root_nodes[i];
							}
						}
						root_nodes.clear();
		maximum_edge_error=0;
		the_surf = 0;
	}
	
	
	/** The Edge2VertexMap is a simple database class.
		The database stores integer values referenced by keys formed by
		a pair of indices. The intended use is to let the key be a pair
		of vertices sharing an edge and the value be a new vertex inserted on
		that edge. */
	class EdgeMap
	{
		typedef pair<int,int>      VertPair;
		typedef map<VertPair, Edge> E2VMap;
		typedef E2VMap::iterator   E2VIter;
		
		E2VMap edge_to_vertex;
		
public:
			
			void set_edge(int ind0, int ind1, const Edge& new_edge)
		{
				edge_to_vertex[VertPair(ind0,ind1)] = new_edge;
		}
		
		bool find_edge(int ind0, int ind1, Edge& edg)
		{
			E2VIter iter = edge_to_vertex.find(VertPair(ind0,ind1));
			
			if(iter == edge_to_vertex.end())
				return false;
			
			edg = iter->second;
			return true;
		}
		
	};
	
	void tessellate(IndexedFaceSet& face_set, ParSurf& s, 
					std::vector<Vec2f> uv_points,
					std::vector<Vec3i> triangles)
	{		
		face_set_ptr = &face_set;
		int i;
		the_surf = &s;
		size_t N = uv_points.size();
		vector<Point> points(N);
		for(i=0;i<N;++i)
		{
			points[i] = create_point(uv_points[i][0], uv_points[i][1]);
			points[i].create_vertex();
		}
		EdgeMap edges;
		vector<Triangle> tris;
		for(size_t j=0;j<triangles.size();++j)
		{
			Edge e[3];
			int v[3];
			
			v[0] = triangles[j][0];
			v[1] = triangles[j][1];
			v[2] = triangles[j][2];
			
			for(int i=0;i<3;++i)
			{
				int k = (i+1)%3;
				if(!edges.find_edge(v[i],v[k],e[i]))
				{
					e[i] = create_edge(points[v[i]], points[v[k]]);
					edges.set_edge(v[k],v[i],e[i].get_opp_edge());
				}
			}
			split(Triangle(e[0],e[1],e[2]));
			
		}
		for(i=0;i<static_cast<int>(root_nodes.size());++i)
		{
			if(root_nodes[i])
			{
				delete root_nodes[i];
			}
		}
		root_nodes.clear();
		maximum_edge_error=0;
		the_surf = 0;
	}
	
}
 

Rev 373	Rev 595
-		1	`/* ----------------------------------------------------------------------- *`
-		2	`* This file is part of GEL, http://www.imm.dtu.dk/GEL`
-		3	`* Copyright (C) the authors and DTU Informatics`
-		4	`* For license and list of authors, see ../../doc/intro.pdf`
-		5	`* ----------------------------------------------------------------------- */`
-		6
1	`#include <map>`	7	`#include <map>`
2	`#include <algorithm>`	8	`#include <algorithm>`
3	`#include <queue>`	9	`#include <queue>`
4		10
5	`#include "Util/Grid2D.h"`	11	`#include "Util/Grid2D.h"`
6	`#include "CGLA/Vec3f.h"`	12	`#include "CGLA/Vec3f.h"`
7	`#include "CGLA/Vec2f.h"`	13	`#include "CGLA/Vec2f.h"`
8		14
9	`#include "tessellate.h"`	15	`#include "tessellate.h"`
10		16
11	`using namespace std;`	17	`using namespace std;`
12	`using namespace CGLA;`	18	`using namespace CGLA;`
13	`using namespace Util;`	19	`using namespace Util;`
14		20
15	`namespace Geometry`	21	`namespace Geometry`
16	`{`	22	`{`
17		23
18	`float MAX_ERR = 5;`	24	`float MAX_ERR = 5;`
19	`float MAX_DIST = 5;`	25	`float MAX_DIST = 5;`
20	`int ADAPTIVE = false;`	26	`int ADAPTIVE = false;`
21		27
22	`#define mymax(x,y) ((x)>(y)?(x):(y))`	28	`#define mymax(x,y) ((x)>(y)?(x):(y))`
23		29
24	`// Test surfaces ---`	30	`// Test surfaces ---`
25		31
26		32
27	`// ----------- CONSTANTS AND GLOBALS ------------`	33	`// ----------- CONSTANTS AND GLOBALS ------------`
28		34
29	`namespace`	35	`namespace`
30	`{`	36	`{`
31	`float maximum_edge_error=0;`	37	`float maximum_edge_error=0;`
32		38
33	`IndexedFaceSet* face_set_ptr = 0;`	39	`IndexedFaceSet* face_set_ptr = 0;`
34	`ParSurf* the_surf = 0;`	40	`ParSurf* the_surf = 0;`
35		41
36	`class TreeNode;`	42	`class TreeNode;`
37	`typedef TreeNode* TreeNodePtr;`	43	`typedef TreeNode* TreeNodePtr;`
38		44
39	`vector<TreeNode*> root_nodes;`	45	`vector<TreeNode*> root_nodes;`
40		46
41	`class Point`	47	`class Point`
42	`{`	48	`{`
43	`friend class Edge;`	49	`friend class Edge;`
44	`friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);`	50	`friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);`
45	`friend float compute_mid_pt_err(const Point&, const Point&, Point&);`	51	`friend float compute_mid_pt_err(const Point&, const Point&, Point&);`
46		52
47		53
48	`Vec2f ppos;`	54	`Vec2f ppos;`
49	`Vec3f pos;`	55	`Vec3f pos;`
50	`mutable int vertex_idx;`	56	`mutable int vertex_idx;`
51		57
52	`public:`	58	`public:`
53		59
54	`Point(const Vec2f& _ppos, const Vec3f& _pos):`	60	`Point(const Vec2f& _ppos, const Vec3f& _pos):`
55	`ppos(_ppos), pos(_pos), vertex_idx(-1)`	61	`ppos(_ppos), pos(_pos), vertex_idx(-1)`
56	`{`	62	`{`
57	`}`	63	`}`
58		64
59	`Point() {}`	65	`Point() {}`
60		66
61	`void create_vertex() const`	67	`void create_vertex() const`
62	`{`	68	`{`
63	`if(vertex_idx == -1)`	69	`if(vertex_idx == -1)`
64	`vertex_idx = face_set_ptr->add_vertex(pos);`	70	`vertex_idx = face_set_ptr->add_vertex(pos);`
65	`}`	71	`}`
66	`};`	72	`};`
67		73
68		74
69		75
70	`inline const Point create_point(float u, float v)`	76	`inline const Point create_point(float u, float v)`
71	`{`	77	`{`
72	`Vec2f ppos(u,v);`	78	`Vec2f ppos(u,v);`
73	`Vec3f pos = (*the_surf)(u,v);`	79	`Vec3f pos = (*the_surf)(u,v);`
74	`return Point(ppos, pos);`	80	`return Point(ppos, pos);`
75	`}`	81	`}`
76		82
77	`float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node);`	83	`float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node);`
78		84
79	`class TreeNode`	85	`class TreeNode`
80	`{`	86	`{`
81	`friend class Edge;`	87	`friend class Edge;`
82	`friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);`	88	`friend float subdiv(int level, const Point&, const Point&, TreeNodePtr&);`
83		89
84	`const Point p;`	90	`const Point p;`
85	`const float error;`	91	`const float error;`
86		92
87	`const TreeNodePtr left;`	93	`const TreeNodePtr left;`
88	`const TreeNodePtr right;`	94	`const TreeNodePtr right;`
89		95
90	`public:`	96	`public:`
91		97
92	`TreeNode(const Point& _p, float _error,`	98	`TreeNode(const Point& _p, float _error,`
93	`TreeNodePtr _left, TreeNodePtr _right):`	99	`TreeNodePtr _left, TreeNodePtr _right):`
94	`p(_p), error(_error), left(_left), right(_right) {}`	100	`p(_p), error(_error), left(_left), right(_right) {}`
95		101
96	`~TreeNode()`	102	`~TreeNode()`
97	`{`	103	`{`
98	`if(left)`	104	`if(left)`
99	`delete left;`	105	`delete left;`
100		106
101	`if(right)`	107	`if(right)`
102	`delete right;`	108	`delete right;`
103	`}`	109	`}`
104		110
105	`};`	111	`};`
106		112
107	`class Edge`	113	`class Edge`
108	`{`	114	`{`
109	`friend Edge create_edge(const Point&, const Point&, bool);`	115	`friend Edge create_edge(const Point&, const Point&, bool);`
110		116
111	`Edge(const Point& _left_point,`	117	`Edge(const Point& _left_point,`
112	`const Point& _right_point,`	118	`const Point& _right_point,`
113	`TreeNodePtr _root,`	119	`TreeNodePtr _root,`
114	`bool _cw):`	120	`bool _cw):`
115	`left_point(_left_point), right_point(_right_point), root(_root), cw(_cw)`	121	`left_point(_left_point), right_point(_right_point), root(_root), cw(_cw)`
116	`{`	122	`{`
117	`}`	123	`}`
118		124
119		125
120	`Point left_point;`	126	`Point left_point;`
121	`Point right_point;`	127	`Point right_point;`
122	`TreeNodePtr root;`	128	`TreeNodePtr root;`
123	`bool cw;`	129	`bool cw;`
124		130
125	`public:`	131	`public:`
126		132
127	`Edge() {}`	133	`Edge() {}`
128		134
129	`float get_error() const`	135	`float get_error() const`
130	`{`	136	`{`
131	`if(root)`	137	`if(root)`
132	`return root->error;`	138	`return root->error;`
133	`return 0.0f;`	139	`return 0.0f;`
134	`}`	140	`}`
135		141
136	`float get_length() const`	142	`float get_length() const`
137	`{`	143	`{`
138	`return (left_point.pos-right_point.pos).length();`	144	`return (left_point.pos-right_point.pos).length();`
139	`}`	145	`}`
140		146
141	`bool is_simple() const`	147	`bool is_simple() const`
142	`{`	148	`{`
143	`if(root->left==0)`	149	`if(root->left==0)`
144	`{`	150	`{`
145	`assert(root->left==0);`	151	`assert(root->left==0);`
146	`assert(root->right==0);`	152	`assert(root->right==0);`
147	`return true;`	153	`return true;`
148	`}`	154	`}`
149	`return false;`	155	`return false;`
150	`}`	156	`}`
151		157
152	`const Edge get_sub_edge(int edge_no) const`	158	`const Edge get_sub_edge(int edge_no) const`
153	`{`	159	`{`
154	`assert(edge_no == 0 \|\| edge_no == 1);`	160	`assert(edge_no == 0 \|\| edge_no == 1);`
155	`assert(!is_simple());`	161	`assert(!is_simple());`
156		162
157	`if(!cw) edge_no = 1-edge_no;`	163	`if(!cw) edge_no = 1-edge_no;`
158		164
159	`if(edge_no ==0)`	165	`if(edge_no ==0)`
160	`return Edge(left_point, root->p, root->left, cw);`	166	`return Edge(left_point, root->p, root->left, cw);`
161	`else`	167	`else`
162	`return Edge(root->p, right_point, root->right, cw);`	168	`return Edge(root->p, right_point, root->right, cw);`
163	`}`	169	`}`
164		170
165	`const Edge get_opp_edge() const`	171	`const Edge get_opp_edge() const`
166	`{`	172	`{`
167	`return Edge(left_point, right_point, root, !cw);`	173	`return Edge(left_point, right_point, root, !cw);`
168	`}`	174	`}`
169		175
170		176
171	`const Point& get_split_point() const`	177	`const Point& get_split_point() const`
172	`{`	178	`{`
173	`assert(!is_simple());`	179	`assert(!is_simple());`
174	`return root->p;`	180	`return root->p;`
175	`}`	181	`}`
176		182
177	`const Point& get_point(int i) const`	183	`const Point& get_point(int i) const`
178	`{`	184	`{`
179	`assert(i==0 \|\| i==1);`	185	`assert(i==0 \|\| i==1);`
180	`if(!cw) i=1-i;`	186	`if(!cw) i=1-i;`
181	`if(i == 0)`	187	`if(i == 0)`
182	`return left_point;`	188	`return left_point;`
183	`return right_point;`	189	`return right_point;`
184	`}`	190	`}`
185		191
186	`int get_point_idx(int i) const`	192	`int get_point_idx(int i) const`
187	`{`	193	`{`
188	`assert(i==0 \|\| i==1);`	194	`assert(i==0 \|\| i==1);`
189	`if(!cw) i=1-i;`	195	`if(!cw) i=1-i;`
190	`if(i == 0)`	196	`if(i == 0)`
191	`return left_point.vertex_idx;`	197	`return left_point.vertex_idx;`
192	`return right_point.vertex_idx;`	198	`return right_point.vertex_idx;`
193	`}`	199	`}`
194		200
195	`};`	201	`};`
196		202
197		203
198	`float compute_mid_pt_err(const Point&, const Point&, Point&);`	204	`float compute_mid_pt_err(const Point&, const Point&, Point&);`
199		205
200	`Edge create_edge(const Point& left_point, const Point& right_point,`	206	`Edge create_edge(const Point& left_point, const Point& right_point,`
201	`bool cw=true)`	207	`bool cw=true)`
202	`{`	208	`{`
203	`TreeNodePtr root=0;`	209	`TreeNodePtr root=0;`
204	`if(ADAPTIVE)`	210	`if(ADAPTIVE)`
205	`subdiv(0, left_point,right_point,root);`	211	`subdiv(0, left_point,right_point,root);`
206	`else`	212	`else`
207	`{`	213	`{`
208	`Point mid_point;`	214	`Point mid_point;`
209	`float err = compute_mid_pt_err(left_point, right_point, mid_point);`	215	`float err = compute_mid_pt_err(left_point, right_point, mid_point);`
210	`root = new TreeNode(mid_point,err,0,0);`	216	`root = new TreeNode(mid_point,err,0,0);`
211	`}`	217	`}`
212	`if(root)`	218	`if(root)`
213	`root_nodes.push_back(root);`	219	`root_nodes.push_back(root);`
214		220
215	`return Edge(left_point, right_point, root, cw);`	221	`return Edge(left_point, right_point, root, cw);`
216	`}`	222	`}`
217		223
218	`class Triangle`	224	`class Triangle`
219	`{`	225	`{`
220	`Edge a,b,c;`	226	`Edge a,b,c;`
221	`public:`	227	`public:`
222	`Triangle(const Edge& _a, const Edge& _b, const Edge& _c):`	228	`Triangle(const Edge& _a, const Edge& _b, const Edge& _c):`
223	`a(_a), b(_b), c(_c) {}`	229	`a(_a), b(_b), c(_c) {}`
224		230
225	`const Edge& operator[](int i) const`	231	`const Edge& operator[](int i) const`
226	`{`	232	`{`
227	`switch(i)`	233	`switch(i)`
228	`{`	234	`{`
229	`case 0: return a;`	235	`case 0: return a;`
230	`case 1: return b;`	236	`case 1: return b;`
231	`case 2: return c;`	237	`case 2: return c;`
232	`}`	238	`}`
233	`return c;`	239	`return c;`
234	`}`	240	`}`
235		241
236	`};`	242	`};`
237		243
238	`float compute_mid_pt_err(const Point& lp, const Point& rp,`	244	`float compute_mid_pt_err(const Point& lp, const Point& rp,`
239	`Point& mid_point)`	245	`Point& mid_point)`
240	`{`	246	`{`
241	`Vec3f m = (lp.pos - rp.pos)/2.0f + rp.pos;`	247	`Vec3f m = (lp.pos - rp.pos)/2.0f + rp.pos;`
242	`Vec2f uvm = (lp.ppos -rp.ppos)/2.0f + rp.ppos;`	248	`Vec2f uvm = (lp.ppos -rp.ppos)/2.0f + rp.ppos;`
243		249
244	`mid_point = create_point(uvm[0], uvm[1]);`	250	`mid_point = create_point(uvm[0], uvm[1]);`
245	`Vec3f diff = mid_point.pos - m;`	251	`Vec3f diff = mid_point.pos - m;`
246	`return diff.length();`	252	`return diff.length();`
247	`}`	253	`}`
248		254
249	`float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node)`	255	`float subdiv(int level, const Point& lp, const Point& rp, TreeNodePtr& node)`
250	`{`	256	`{`
251	`Point mid_point;`	257	`Point mid_point;`
252	`float err = compute_mid_pt_err(lp,rp,mid_point);`	258	`float err = compute_mid_pt_err(lp,rp,mid_point);`
253		259
254	`Vec3f dist = lp.pos-rp.pos;`	260	`// Vec3f dist = lp.pos-rp.pos;`
255	`if (level>10)`	261	`if (level>10)`
256	`{`	262	`{`
257	`node = new TreeNode(mid_point,err,0,0);`	263	`node = new TreeNode(mid_point,err,0,0);`
258	`return err;`	264	`return err;`
259	`}`	265	`}`
260		266
261	`TreeNodePtr lnode;`	267	`TreeNodePtr lnode;`
262	`float lerr = subdiv(level+1,lp, mid_point, lnode);`	268	`float lerr = subdiv(level+1,lp, mid_point, lnode);`
263		269
264	`TreeNodePtr rnode;`	270	`TreeNodePtr rnode;`
265	`float rerr = subdiv(level+1,mid_point, rp, rnode);`	271	`float rerr = subdiv(level+1,mid_point, rp, rnode);`
266		272
267	`float max_err = mymax(err, mymax(lerr,rerr));`	273	`float max_err = mymax(err, mymax(lerr,rerr));`
268		274
269	`if(max_err < MAX_ERR)`	275	`if(max_err < MAX_ERR)`
270	`{`	276	`{`
271	`if(lnode) delete lnode;`	277	`if(lnode) delete lnode;`
272	`if(rnode) delete rnode;`	278	`if(rnode) delete rnode;`
273		279
274	`node = new TreeNode(mid_point, max_err, 0,0);`	280	`node = new TreeNode(mid_point, max_err, 0,0);`
275	`}`	281	`}`
276	`else`	282	`else`
277	`{`	283	`{`
278	`mid_point.create_vertex();`	284	`mid_point.create_vertex();`
279	`node = new TreeNode(mid_point, max_err, lnode,rnode);`	285	`node = new TreeNode(mid_point, max_err, lnode,rnode);`
280	`}`	286	`}`
281		287
282	`return max_err;`	288	`return max_err;`
283	`}`	289	`}`
284		290
285		291
286		292
287	`void split(const Triangle orig_tri)`	293	`void split(const Triangle orig_tri)`
288	`{`	294	`{`
289	`queue<Triangle> tri_queue;`	295	`queue<Triangle> tri_queue;`
290	`tri_queue.push(orig_tri);`	296	`tri_queue.push(orig_tri);`
291		297
292		298
293	`while(!tri_queue.empty())`	299	`while(!tri_queue.empty())`
294	`{`	300	`{`
295	`const Triangle tri = tri_queue.front();`	301	`const Triangle tri = tri_queue.front();`
296	`tri_queue.pop();`	302	`tri_queue.pop();`
297		303
298	`int no_split_edges = 0;`	304	`int no_split_edges = 0;`
299		305
300	`for(int i=0;i<3;++i)`	306	`for(int i=0;i<3;++i)`
301	`if(!tri[i].is_simple())`	307	`if(!tri[i].is_simple())`
302	`++no_split_edges;`	308	`++no_split_edges;`
303		309
304	`if(no_split_edges==0 /* \|\| level > 2000*/)`	310	`if(no_split_edges==0 /* \|\| level > 2000*/)`
305	`{`	311	`{`
306	`int a = tri[0].get_point_idx(0);`	312	`int a = tri[0].get_point_idx(0);`
307	`int b = tri[1].get_point_idx(0);`	313	`int b = tri[1].get_point_idx(0);`
308	`int c = tri[2].get_point_idx(0);`	314	`int c = tri[2].get_point_idx(0);`
309	`face_set_ptr->add_face(Vec3i(a,b,c));`	315	`face_set_ptr->add_face(Vec3i(a,b,c));`
310	`maximum_edge_error = mymax(mymax(tri[0].get_error(),`	316	`maximum_edge_error = mymax(mymax(tri[0].get_error(),`
311	`tri[1].get_error()),`	317	`tri[1].get_error()),`
312	`tri[2].get_error());`	318	`tri[2].get_error());`
313	`}`	319	`}`
314	`else if(no_split_edges==3)`	320	`else if(no_split_edges==3)`
315	`{`	321	`{`
316	`const Point& p0 = tri[0].get_point(0);`	322	`const Point& p0 = tri[0].get_point(0);`
317	`const Point& p1 = tri[1].get_point(0);`	323	`const Point& p1 = tri[1].get_point(0);`
318	`const Point& p2 = tri[2].get_point(0);`	324	`const Point& p2 = tri[2].get_point(0);`
319		325
320	`const Point sp[3] = {tri[0].get_split_point(),`	326	`const Point sp[3] = {tri[0].get_split_point(),`
321	`tri[1].get_split_point(),`	327	`tri[1].get_split_point(),`
322	`tri[2].get_split_point()};`	328	`tri[2].get_split_point()};`
323		329
324	`Edge sp_edg[3] = {create_edge(sp[0],p2),`	330	`Edge sp_edg[3] = {create_edge(sp[0],p2),`
325	`create_edge(sp[1],p0),`	331	`create_edge(sp[1],p0),`
326	`create_edge(sp[2],p1)};`	332	`create_edge(sp[2],p1)};`
327		333
328	`int i_min=0;`	334	`int i_min=0;`
329	`float l_min = sp_edg[0].get_length();`	335	`float l_min = sp_edg[0].get_length();`
330	`int i;`	336	`int i;`
331	`for(i=1;i<3;++i)`	337	`for(i=1;i<3;++i)`
332	`{`	338	`{`
333	`float len = sp_edg[i].get_length();`	339	`float len = sp_edg[i].get_length();`
334	`if(len<l_min)`	340	`if(len<l_min)`
335	`{`	341	`{`
336	`l_min = len;`	342	`l_min = len;`
337	`i_min = i;`	343	`i_min = i;`
338	`}`	344	`}`
339	`}`	345	`}`
340		346
341	`i = i_min;`	347	`i = i_min;`
342	`int j = (i+1)%3;`	348	`int j = (i+1)%3;`
343	`int k = (i+2)%3;`	349	`int k = (i+2)%3;`
344		350
345	`Edge spi_spj = create_edge(sp[i],sp[j]);`	351	`Edge spi_spj = create_edge(sp[i],sp[j]);`
346	`Edge spi_spk = create_edge(sp[i],sp[k]);`	352	`Edge spi_spk = create_edge(sp[i],sp[k]);`
347		353
348	`Triangle tri0(tri[i].get_sub_edge(1),`	354	`Triangle tri0(tri[i].get_sub_edge(1),`
349	`tri[j].get_sub_edge(0),`	355	`tri[j].get_sub_edge(0),`
350	`spi_spj.get_opp_edge());`	356	`spi_spj.get_opp_edge());`
351	`Triangle tri1(spi_spj,`	357	`Triangle tri1(spi_spj,`
352	`tri[j].get_sub_edge(1),`	358	`tri[j].get_sub_edge(1),`
353	`sp_edg[i].get_opp_edge());`	359	`sp_edg[i].get_opp_edge());`
354	`Triangle tri2(sp_edg[i],`	360	`Triangle tri2(sp_edg[i],`
355	`tri[k].get_sub_edge(0),`	361	`tri[k].get_sub_edge(0),`
356	`spi_spk.get_opp_edge());`	362	`spi_spk.get_opp_edge());`
357	`Triangle tri3(spi_spk,`	363	`Triangle tri3(spi_spk,`
358	`tri[k].get_sub_edge(1),`	364	`tri[k].get_sub_edge(1),`
359	`tri[i].get_sub_edge(0));`	365	`tri[i].get_sub_edge(0));`
360		366
361	`tri_queue.push(tri0);`	367	`tri_queue.push(tri0);`
362	`tri_queue.push(tri1);`	368	`tri_queue.push(tri1);`
363	`tri_queue.push(tri2);`	369	`tri_queue.push(tri2);`
364	`tri_queue.push(tri3);`	370	`tri_queue.push(tri3);`
365	`}`	371	`}`
366	`else if(no_split_edges==2)`	372	`else if(no_split_edges==2)`
367	`{`	373	`{`
368	`int i;`	374	`int i;`
369	`if(tri[2].is_simple()) i = 0;`	375	`if(tri[2].is_simple()) i = 0;`
370	`else if(tri[1].is_simple()) i=2;`	376	`else if(tri[1].is_simple()) i=2;`
371	`else i = 1;`	377	`else i = 1;`
372		378
373	`int j = (i+1)%3;`	379	`int j = (i+1)%3;`
374	`int k = (i+2)%3;`	380	`int k = (i+2)%3;`
375		381
376	`const Point& spi = tri[i].get_split_point();`	382	`const Point& spi = tri[i].get_split_point();`
377	`const Point& spj = tri[j].get_split_point();`	383	`const Point& spj = tri[j].get_split_point();`
378		384
379	`const Point& pi = tri[i].get_point(0);`	385	`const Point& pi = tri[i].get_point(0);`
380	`const Point& pk = tri[k].get_point(0);`	386	`const Point& pk = tri[k].get_point(0);`
381		387
382	`Edge spi_spj = create_edge(spi, spj);`	388	`Edge spi_spj = create_edge(spi, spj);`
383	`Edge pi_spj = create_edge(pi, spj);`	389	`Edge pi_spj = create_edge(pi, spj);`
384	`Edge pk_spi = create_edge(pk, spi);`	390	`Edge pk_spi = create_edge(pk, spi);`
385		391
386	`Triangle tri0(tri[i].get_sub_edge(1),`	392	`Triangle tri0(tri[i].get_sub_edge(1),`
387	`tri[j].get_sub_edge(0),`	393	`tri[j].get_sub_edge(0),`
388	`spi_spj.get_opp_edge());`	394	`spi_spj.get_opp_edge());`
389		395
390		396
391	`if(pi_spj.get_length() < pk_spi.get_length())`	397	`if(pi_spj.get_length() < pk_spi.get_length())`
392	`{`	398	`{`
393	`Triangle tri1(spi_spj,`	399	`Triangle tri1(spi_spj,`
394	`pi_spj.get_opp_edge(),`	400	`pi_spj.get_opp_edge(),`
395	`tri[i].get_sub_edge(0));`	401	`tri[i].get_sub_edge(0));`
396	`Triangle tri2(pi_spj,`	402	`Triangle tri2(pi_spj,`
397	`tri[j].get_sub_edge(1),`	403	`tri[j].get_sub_edge(1),`
398	`tri[k]);`	404	`tri[k]);`
399	`tri_queue.push(tri1);`	405	`tri_queue.push(tri1);`
400	`tri_queue.push(tri2);`	406	`tri_queue.push(tri2);`
401	`}`	407	`}`
402	`else`	408	`else`
403	`{`	409	`{`
404	`Triangle tri1(tri[i].get_sub_edge(0),`	410	`Triangle tri1(tri[i].get_sub_edge(0),`
405	`pk_spi.get_opp_edge(),`	411	`pk_spi.get_opp_edge(),`
406	`tri[k]);`	412	`tri[k]);`
407	`Triangle tri2(spi_spj,`	413	`Triangle tri2(spi_spj,`
408	`tri[j].get_sub_edge(1),`	414	`tri[j].get_sub_edge(1),`
409	`pk_spi);`	415	`pk_spi);`
410	`tri_queue.push(tri1);`	416	`tri_queue.push(tri1);`
411	`tri_queue.push(tri2);`	417	`tri_queue.push(tri2);`
412	`}`	418	`}`
413	`tri_queue.push(tri0);`	419	`tri_queue.push(tri0);`
414	`}`	420	`}`
415	`else if(no_split_edges==1)`	421	`else if(no_split_edges==1)`
416	`{`	422	`{`
417	`int i;`	423	`int i;`
418	`if(!tri[0].is_simple()) i=0;`	424	`if(!tri[0].is_simple()) i=0;`
419	`else if(!tri[1].is_simple()) i=1;`	425	`else if(!tri[1].is_simple()) i=1;`
420	`else i= 2;`	426	`else i= 2;`
421		427
422	`int j = (i+1)%3;`	428	`int j = (i+1)%3;`
423	`int k = (i+2)%3;`	429	`int k = (i+2)%3;`
424		430
425	`const Point& pk = tri[k].get_point(0);`	431	`const Point& pk = tri[k].get_point(0);`
426	`const Point& spi = tri[i].get_split_point();`	432	`const Point& spi = tri[i].get_split_point();`
427		433
428	`Edge spi_pk = create_edge(spi, pk);`	434	`Edge spi_pk = create_edge(spi, pk);`
429		435
430	`Triangle tri0(tri[i].get_sub_edge(1),`	436	`Triangle tri0(tri[i].get_sub_edge(1),`
431	`tri[j],`	437	`tri[j],`
432	`spi_pk.get_opp_edge());`	438	`spi_pk.get_opp_edge());`
433		439
434	`Triangle tri1(tri[i].get_sub_edge(0),`	440	`Triangle tri1(tri[i].get_sub_edge(0),`
435	`spi_pk,`	441	`spi_pk,`
436	`tri[k]);`	442	`tri[k]);`
437		443
438	`tri_queue.push(tri0);`	444	`tri_queue.push(tri0);`
439	`tri_queue.push(tri1);`	445	`tri_queue.push(tri1);`
440		446
441	`}`	447	`}`
442	`}`	448	`}`
443	`}`	449	`}`
444	`}`	450	`}`
445		451
446		452
447	`void tessellate(IndexedFaceSet& face_set, ParSurf& s, Grid2D<Vec3f>& inigrid)`	453	`void tessellate(IndexedFaceSet& face_set, ParSurf& s, Grid2D<Vec3f>& inigrid)`
448	`{`	454	`{`
449	`face_set_ptr = &face_set;`	455	`face_set_ptr = &face_set;`
450	`int i;`	456	`int i;`
451	`int n=inigrid.get_xdim()-1;`	457	`int n=inigrid.get_xdim()-1;`
452	`int m=inigrid.get_ydim()-1;`	458	`int m=inigrid.get_ydim()-1;`
453	`the_surf = &s;`	459	`the_surf = &s;`
454	`Grid2D<Point> points(n+1,m+1);`	460	`Grid2D<Point> points(n+1,m+1);`
455	`for(i=0;i<=n;++i)`	461	`for(i=0;i<=n;++i)`
456	`for(int j=0;j<=m;++j)`	462	`for(int j=0;j<=m;++j)`
457	`{`	463	`{`
458	`Vec3f p = inigrid(i,j);`	464	`Vec3f p = inigrid(i,j);`
459	`points(i,j) = create_point(p[0], p[1]);`	465	`points(i,j) = create_point(p[0], p[1]);`
460	`points(i,j).create_vertex();`	466	`points(i,j).create_vertex();`
461	`}`	467	`}`
462		468
463		469
464	`Grid2D<Edge> h_edges(n,m+1);`	470	`Grid2D<Edge> h_edges(n,m+1);`
465	`Grid2D<Edge> v_edges(n+1,m);`	471	`Grid2D<Edge> v_edges(n+1,m);`
466	`for(i=0;i<=n;++i)`	472	`for(i=0;i<=n;++i)`
467	`for(int j=0;j<=m;++j)`	473	`for(int j=0;j<=m;++j)`
468	`{`	474	`{`
469	`if(i<n)`	475	`if(i<n)`
470	`h_edges(i,j) = create_edge(points(i,j), points(i+1,j));`	476	`h_edges(i,j) = create_edge(points(i,j), points(i+1,j));`
471	`if(j<m)`	477	`if(j<m)`
472	`v_edges(i,j) = create_edge(points(i,j), points(i,j+1));`	478	`v_edges(i,j) = create_edge(points(i,j), points(i,j+1));`
473	`}`	479	`}`
474		480
475	`for(i=0;i<n;++i)`	481	`for(i=0;i<n;++i)`
476	`for(int j=0;j<m;++j)`	482	`for(int j=0;j<m;++j)`
477	`{`	483	`{`
478	`Edge diag = create_edge(points(i,j), points(i+1,j+1));`	484	`Edge diag = create_edge(points(i,j), points(i+1,j+1));`
479	`Triangle tri0(h_edges(i,j), v_edges(i+1,j), diag.get_opp_edge());`	485	`Triangle tri0(h_edges(i,j), v_edges(i+1,j), diag.get_opp_edge());`
480	`Triangle tri1(diag, h_edges(i,j+1).get_opp_edge(),`	486	`Triangle tri1(diag, h_edges(i,j+1).get_opp_edge(),`
481	`v_edges(i,j).get_opp_edge());`	487	`v_edges(i,j).get_opp_edge());`
482	`split(tri0);`	488	`split(tri0);`
483	`split(tri1);`	489	`split(tri1);`
484	`}`	490	`}`
485	`for(i=0;i<static_cast<int>(root_nodes.size());++i)`	491	`for(i=0;i<static_cast<int>(root_nodes.size());++i)`
486	`{`	492	`{`
487	`if(root_nodes[i])`	493	`if(root_nodes[i])`
488	`{`	494	`{`
489	`delete root_nodes[i];`	495	`delete root_nodes[i];`
490	`}`	496	`}`
491	`}`	497	`}`
492	`root_nodes.clear();`	498	`root_nodes.clear();`
493	`maximum_edge_error=0;`	499	`maximum_edge_error=0;`
494	`the_surf = 0;`	500	`the_surf = 0;`

Subversion Repositories gelsvn

(root)/trunk/src/Geometry/tessellate.cpp – Rev 373 → 595