Subversion Repositories gelsvn

Rev

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

Rev Author Line No. Line
2 bj 1
#include "Vec3f.h"
2
#include "Quaternion.h"
3
 
12 jab 4
namespace CGLA 
5
{
171 bj 6
 
7
	void Quaternion::make_rot(const Mat4x4f& m)
8
	{
9
		float trace = m[0][0] + m[1][1] + m[2][2]  + 1.0f;
10
 
11
		//If the trace of the matrix is greater than zero, then
12
		//perform an "instant" calculation.
13
		if ( trace > TINY ) 
14
		{
15
			float S = sqrt(trace) * 2.0f;
16
			qv = Vec3f(m[2][1] - m[1][2], m[0][2] - m[2][0], m[1][0] - m[0][1] );
17
			qv /= S;
18
			qw = 0.25f * S;
19
		}
20
		else
21
		{
22
			//If the trace of the matrix is equal to zero (or negative...) then identify
23
			//which major diagonal element has the greatest value.
24
			//Depending on this, calculate the following:
25
 
26
			if ( m[0][0] > m[1][1] && m[0][0] > m[2][2] )  {	// Column 0: 
27
				float S  = sqrt( 1.0 + m[0][0] - m[1][1] - m[2][2] ) * 2.0f;
28
				qv[0] = 0.25f * S;
29
				qv[1] = (m[1][0] + m[0][1] ) / S;
30
				qv[2] = (m[0][2] + m[2][0] ) / S;
31
				qw = (m[2][1] - m[1][3] ) / S;
32
			} else if ( m[1][1] > m[2][2] ) {			// Column 1: 
33
				float S  = sqrt( 1.0 + m[1][1] - m[0][0] - m[2][2] ) * 2.0f;
34
				qv[0] = (m[1][0] + m[0][1] ) / S;
35
				qv[1] = 0.25f * S;
36
				qv[2] = (m[2][1] + m[1][2] ) / S;
37
				qw = (m[0][2] - m[2][0] ) / S;
38
			} else {						// Column 2:
39
				float S  = sqrt( 1.0 + m[2][2] - m[0][0] - m[1][1] ) * 2.0f;
40
				qv[0] = (m[0][2] + m[2][0] ) / S;
41
				qv[1] = (m[2][1] + m[1][2] ) / S;
42
				qv[2] = 0.25f * S;
43
				qw = (m[1][0] - m[0][1] ) / S;
44
			}
45
		}
46
		//The quaternion is then defined as:
47
		//  Q = | X Y Z W |
48
	}
49
 
50
	void Quaternion::make_rot(float angle, const Vec3f& v)
12 jab 51
  {
52
    angle = angle/2;
53
    qv = CGLA::normalize(v);
54
    qv *= sin(angle);
55
    qw  = cos(angle);
56
  };
2 bj 57
 
12 jab 58
  void Quaternion::make_rot(const Vec3f& _v0, const Vec3f& _v1)
59
  {
60
    Vec3f v0 = CGLA::normalize(_v0);
61
    Vec3f v1 = CGLA::normalize(_v1);
62
    qv = cross(v0, v1);
63
    float l = qv.length();
64
    if(l<TINY)
65
      qv = Vec3f(1,0,0);
66
    else
67
      qv.normalize();
68
    float a = acos(dot(v0,v1))/2;
69
    qw  = cos(a);
70
    qv *= sin(a);	
71
  };
2 bj 72
 
12 jab 73
  void Quaternion::get_rot(float& angle, Vec3f& v) 
74
  {
75
    angle=2*acos(qw);
2 bj 76
 
12 jab 77
    if (angle<TINY) 
78
      v = Vec3f(1,0,0);
79
    else 
80
      v = qv / sin(angle);
2 bj 81
 
12 jab 82
    if (angle>M_PI)
83
      v = -v;
2 bj 84
 
12 jab 85
    v.normalize();
86
  }
2 bj 87
 
12 jab 88
  Mat3x3f Quaternion::get_mat3x3f() const
89
  {
90
    float s=2/norm();
91
    float m[9] = {1 - s*(qv[1]*qv[1]+qv[2]*qv[2]), 
92
		  s*(qv[0]*qv[1]-qw*qv[2]), 
93
		  s*(qv[0]*qv[2]+qw*qv[1]), 
94
		  s*(qv[0]*qv[1]+qw*qv[2]), 
95
		  1 - s*(qv[0]*qv[0]+qv[2]*qv[2]), 
96
		  s*(qv[1]*qv[2]-qw*qv[0]), 
97
		  s*(qv[0]*qv[2]-qw*qv[1]), 
98
		  s*(qv[1]*qv[2]+qw*qv[0]), 
99
		  1 - s*(qv[0]*qv[0]+qv[1]*qv[1])};
100
    Mat3x3f mat;
101
    raw_assign(mat, m);
102
    return mat;
103
  }
2 bj 104
 
105
 
12 jab 106
  //This function just need to call the right initialiser
2 bj 107
 
12 jab 108
  Mat4x4f Quaternion::get_mat4x4f() const
109
  {
110
    float s=2/norm();
111
    float m[16] = {1 - s*(qv[1]*qv[1]+qv[2]*qv[2]), 
112
		   s*(qv[0]*qv[1]-qw*qv[2]), 
113
		   s*(qv[0]*qv[2]+qw*qv[1]), 
114
		   float(0),
115
		   s*(qv[0]*qv[1]+qw*qv[2]), 
116
		   1 - s*(qv[0]*qv[0]+qv[2]*qv[2]), 
117
		   s*(qv[1]*qv[2]-qw*qv[0]), 
118
		   float(0),
119
		   s*(qv[0]*qv[2]-qw*qv[1]), 
120
		   s*(qv[1]*qv[2]+qw*qv[0]), 
121
		   1 - s*(qv[0]*qv[0]+qv[1]*qv[1]), 
122
		   float(0),
123
		   float(0),             
124
		   float(0),                   
125
		   float(0),                
126
		   float(1)};
127
    Mat4x4f mat;
128
    raw_assign(mat, m);
129
    return mat;
130
  }
2 bj 131
 
132
 
12 jab 133
  Vec3f Quaternion::apply(const Vec3f& vec) const
134
  {
135
    return Vec3f((*this)*Quaternion(vec)*inverse());
136
  }
2 bj 137
 
138
 
139
}