Subversion Repositories gelsvn

Rev

Rev 214 | Rev 315 | Go to most recent revision | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
39 bj 1
#ifndef __QEM_H
2
#define __QEM_H
3
 
4
#include <cfloat>
5
#include "CGLA/Vec3d.h"
6
#include "CGLA/Mat3x3d.h"
7
 
89 jab 8
 
39 bj 9
namespace
10
{
11
	inline const CGLA::Mat3x3d direct_product(const CGLA::Vec3d& v0, const CGLA::Vec3d& v1)
12
		{
13
			CGLA::Mat3x3d m;
14
			for(int i=0;i<3;++i)
15
				for(int j=0;j<3;++j)
16
					m[i][j] = v0[i]*v1[j];
17
			return m;
18
		}
19
}
20
 
21
namespace Geometry
22
{
23
	class QEM
24
		{
25
			CGLA::Mat3x3d A;
26
			CGLA::Vec3d   b;
27
			double   c;
28
		public:
29
 
136 jab 30
			QEM(): A(0), b(0), c(0) {}
39 bj 31
 
32
			QEM(const CGLA::Vec3d& p0, const CGLA::Vec3d& n0, double w=1.0f):
33
				A(direct_product(n0,n0) * w), 
34
				b(-n0*dot(n0,p0) * w), 
35
				c(dot(p0,n0)*dot(p0,n0) * w) {}
36
 
37
 
38
			void operator+=(const QEM& q)
39
				{
40
					A += q.A;
41
					b += q.b;
42
					c += q.c;
43
				}
44
 
220 jab 45
			float error(const CGLA::Vec3d& p) const
39 bj 46
				{
47
					return dot(p,A*p) + 2*dot(b,p)+ c;
48
				}
49
 
50
			double determinant() const
51
				{
52
					return CGLA::determinant(A);
53
				}
54
 
214 jab 55
			const CGLA::Vec3d grad(const CGLA::Vec3d& _p) const
39 bj 56
				{
57
					CGLA::Vec3d p(_p);
214 jab 58
					CGLA::Vec3d g(2*A*p+2*b);
39 bj 59
					return g;
60
				}
61
 
214 jab 62
			CGLA::Vec3d opt_pos(double QEM_thresh = 0.005) const;
39 bj 63
 
64
		};
65
}
66
 
67
namespace GEO = Geometry;
68
 
69
#endif