manolas
/
vcglib
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

* change all remaining Transpose to transpose, 

* update the gl/math wrappers to make them more Eigen friendly
  (and remove the useless, and not used, and somehow dangerous
   *Direct and *E functions)
* add automatic reinterpret_casting from Eigen::Matrix to vcg
  specialized types
This commit is contained in:
Paolo Cignoni 2008-10-29 00:05:44 +00:00
parent 0361427bc0
commit 07f2e976ea
15 changed files with 533 additions and 377 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

36
apps/ptx2ply/ptx2ply.cpp
View File

@ -43,7 +43,7 @@ using namespace tri;
class MyFace; class MyFace;
class MyEdgeC; class MyEdgeC;
class MyFaceC; class MyFaceC;
class MyVertexC:public VertexVCVN<float,MyEdgeC,MyFace>{}; class MyVertexC:public VertexVCVN<float,MyEdgeC,MyFace>{};
class MyFaceC :public FaceFN<MyVertexC,MyEdgeC,MyFaceC>{}; class MyFaceC :public FaceFN<MyVertexC,MyEdgeC,MyFaceC>{};
class MyMeshC: public tri::TriMesh< std::vector<MyVertexC>, std::vector<MyFaceC > >{}; class MyMeshC: public tri::TriMesh< std::vector<MyVertexC>, std::vector<MyFaceC > >{};
@ -148,11 +148,11 @@ int readmesh(FILE* fp)
if(hascolor && savecolor) if(hascolor && savecolor)
{ {
viC = Allocator<MyMeshC>::AddVertices(currentmeshC,(rownum*colnum)); viC = Allocator<MyMeshC>::AddVertices(currentmeshC,(rownum*colnum));
} }
else else
{ {
vi = Allocator<MyMesh>::AddVertices(currentmesh,(rownum*colnum)); vi = Allocator<MyMesh>::AddVertices(currentmesh,(rownum*colnum));
} }
// parse the first line.... // parse the first line....
@ -360,7 +360,7 @@ int readmesh(FILE* fp)
if(hascolor && savecolor) if(hascolor && savecolor)
printf("\nV: %8i F: %8i \n", currentmeshC.vn, currentmeshC.fn); printf("\nV: %8i F: %8i \n", currentmeshC.vn, currentmeshC.fn);
else else
printf("\nV: %8i F: %8i \n", currentmesh.vn, currentmesh.fn); printf("\nV: %8i F: %8i \n", currentmesh.vn, currentmesh.fn);
} }
@ -419,7 +419,7 @@ int readmesh(FILE* fp)
if(hascolor && savecolor) if(hascolor && savecolor)
printf("V: %8i F: %8i \n", currentmeshC.vn, currentmeshC.fn); printf("V: %8i F: %8i \n", currentmeshC.vn, currentmeshC.fn);
else else
printf("V: %8i F: %8i \n", currentmesh.vn, currentmesh.fn); printf("V: %8i F: %8i \n", currentmesh.vn, currentmesh.fn);
@ -464,7 +464,7 @@ int readmesh(FILE* fp)
} }
} }
Transpose(currtrasf); currtrasf = currtrasf.transpose().eval();
// apply tranformation // apply tranformation
if(hascolor && savecolor) if(hascolor && savecolor)
@ -496,7 +496,7 @@ int readmesh(FILE* fp)
if(! onlypoints) if(! onlypoints)
int unref = tri::Clean<MyMesh>::RemoveUnreferencedVertex(currentmesh); int unref = tri::Clean<MyMesh>::RemoveUnreferencedVertex(currentmesh);
} }
return 0; return 0;
} }
@ -568,7 +568,7 @@ void dounpack(FILE* fp)
} }
// skip a mesh // skip a mesh
int skipmesh(FILE* fp) int skipmesh(FILE* fp)
{ {
int colnum; int colnum;
@ -595,7 +595,7 @@ int skipmesh(FILE* fp)
fread(&linebuf,1,1,fp); fread(&linebuf,1,1,fp);
while(linebuf != '\n') while(linebuf != '\n')
fread(&linebuf,1,1,fp); fread(&linebuf,1,1,fp);
} }
return 0; return 0;
} }
@ -643,17 +643,17 @@ void parseparams(int argn, char** argvect)
dumpit = true; dumpit = true;
printf("dumping # %i chars\n",todump); printf("dumping # %i chars\n",todump);
} }
if(argvect[pit][1] == 'u') // unpack the file in different if(argvect[pit][1] == 'u') // unpack the file in different
{ {
unpack = true; unpack = true;
printf("UNPACKING \n"); printf("UNPACKING \n");
} }
if(argvect[pit][1] == 'c') // save color if present if(argvect[pit][1] == 'c') // save color if present
{ {
savecolor = true; savecolor = true;
printf("saving color \n"); printf("saving color \n");
} }
if(argvect[pit][1] == 'k') // keep all if(argvect[pit][1] == 'k') // keep all
{ {
saveall = true; saveall = true;
printf("keeping all elements \n"); printf("keeping all elements \n");
@ -673,7 +673,7 @@ void parseparams(int argn, char** argvect)
switchside = true; switchside = true;
printf("swapped triangulation \n"); printf("swapped triangulation \n");
} }
} }
} }
@ -750,19 +750,19 @@ int main(int argc, char *argv[])
dumpit = false; dumpit = false;
unpack = false; unpack = false;
savecolor = false; savecolor = false;
saveall = false; saveall = false;
flipfaces = false; flipfaces = false;
onlypoints = false; onlypoints = false;
switchside = false; switchside = false;
if(argc < 2) if(argc < 2)
printhelp(); printhelp();
//-- //--
parseparams(argc, argv); parseparams(argc, argv);
strcpy(modelname,argv[1]); strcpy(modelname,argv[1]);
modelname[strlen(argv[1])-4] = '\0'; modelname[strlen(argv[1])-4] = '\0';
@ -794,7 +794,7 @@ int main(int argc, char *argv[])
fclose(fp); fclose(fp);
exit(0); exit(0);
} }
printf("reading "); printf("reading ");
readmesh(fp); readmesh(fp);

29
apps/test/quaternion/quat_test.cpp
View File

@ -47,13 +47,13 @@ ostream &operator<<(ostream &o, Matrix44f &m) {
return o; return o;
} }
bool verify(Quaternionf &q){ bool verify(Quaternionf &q){
cout << "Quaternion: " << q << endl; cout << "Quaternion: " << q << endl;
Matrix33f m; Matrix33f m;
q.ToMatrix(m); q.ToMatrix(m);
cout << "To Matrix: " << m << endl; cout << "To Matrix: " << m << endl;
cout << "Row norms: " << m.GetRow(0).Norm() << " " cout << "Row norms: " << m.GetRow(0).Norm() << " "
<< m.GetRow(1).Norm() << " " << m.GetRow(1).Norm() << " "
<< m.GetRow(2).Norm() << endl; << m.GetRow(2).Norm() << endl;
Point3f p(3, 4, 5); Point3f p(3, 4, 5);
Point3f qp = q.Rotate(p); Point3f qp = q.Rotate(p);
@ -68,14 +68,13 @@ bool verify(Quaternionf &q){
bool verify(Matrix33f &m) { bool verify(Matrix33f &m) {
cout << "Matrix: " << m << endl; cout << "Matrix: " << m << endl;
cout << "Det: " << m.Determinant() << endl; cout << "Det: " << m.Determinant() << endl;
cout << "Row norms: " << m.GetRow(0).Norm() << " " cout << "Row norms: " << m.GetRow(0).Norm() << " "
<< m.GetRow(1).Norm() << " " << m.GetRow(1).Norm() << " "
<< m.GetRow(2).Norm() << endl; << m.GetRow(2).Norm() << endl;
cout << "Column norms: " << m.GetColumn(0).Norm() << " " cout << "Column norms: " << m.GetColumn(0).Norm() << " "
<< m.GetColumn(1).Norm() << " " << m.GetColumn(1).Norm() << " "
<< m.GetColumn(2).Norm() << endl; << m.GetColumn(2).Norm() << endl;
Matrix33f im = m; Matrix33f im = m.transpose() * m;
im.Transpose();
im = im*m; im = im*m;
cout << "Check ortonormality: " << im << endl; cout << "Check ortonormality: " << im << endl;
Quaternionf q; Quaternionf q;
@ -87,7 +86,7 @@ bool verify(Matrix33f &m) {
cout << "Norm: " << axis.SquaredNorm() + q[0]*q[0] << endl; cout << "Norm: " << axis.SquaredNorm() + q[0]*q[0] << endl;
axis.Normalize(); axis.Normalize();
cout << "angle: " << 2*acos(q[0]) << " Axis: " << axis << endl; cout << "angle: " << 2*acos(q[0]) << " Axis: " << axis << endl;
Point3f p(3, 4, 5); Point3f p(3, 4, 5);
Point3f qp = q.Rotate(p); Point3f qp = q.Rotate(p);
Point3f mp = m*p; Point3f mp = m*p;
@ -112,8 +111,8 @@ int main() {
verify(q); verify(q);
Matrix33f m; Matrix33f m;
m[0][0] = 0.70145; m[0][1] = 0.372035; m[0][2] = 0.607913; m[0][0] = 0.70145; m[0][1] = 0.372035; m[0][2] = 0.607913;
m[1][0] = -0.628023; m[1][1] = 0.725922; m[1][2] = 0.2804; m[1][0] = -0.628023; m[1][1] = 0.725922; m[1][2] = 0.2804;
m[2][0] = -0.336978; m[2][1] = -0.57847; m[2][2] = 0.742845; m[2][0] = -0.336978; m[2][1] = -0.57847; m[2][2] = 0.742845;
cout << "verify matrix: " << endl; cout << "verify matrix: " << endl;
@ -129,7 +128,7 @@ int main() {
//Matrix33f m; //Matrix33f m;
cout << "matrix: " << m << endl; cout << "matrix: " << m << endl;
Point3f p(3, 4, 5); Point3f p(3, 4, 5);
//test this matrix: //test this matrix:
cout << "norms: " << m.GetRow(0).Norm() << " " << m.GetRow(1).Norm() << " " << m.GetRow(2).Norm() << endl; cout << "norms: " << m.GetRow(0).Norm() << " " << m.GetRow(1).Norm() << " " << m.GetRow(2).Norm() << endl;
@ -143,14 +142,14 @@ Point3f p(3, 4, 5);
} }
cout << endl; cout << endl;
} }
cout << endl; cout << endl;
cout << "Point: " << p[0] << " " << p[1] << " " << p[2] << endl; cout << "Point: " << p[0] << " " << p[1] << " " << p[2] << endl;
Point3f r = q.Rotate(p); Point3f r = q.Rotate(p);
cout << "Rotated by q: " << r[0] << " " << r[1] << " " << r[2] << endl; cout << "Rotated by q: " << r[0] << " " << r[1] << " " << r[2] << endl;
r = m*p; r = m*p;
cout << "Rotated by m: " << r[0] << " " << r[1] << " " << r[2] << endl; cout << "Rotated by m: " << r[0] << " " << r[1] << " " << r[2] << endl;
q.FromMatrix(m); q.FromMatrix(m);
cout << "quaternion: " << q[0] << " " << q[1] << " " << q[2] << " " << q[3] <<endl; cout << "quaternion: " << q[0] << " " << q[1] << " " << q[2] << " " << q[3] <<endl;

287
vcg/math/deprecated_matrix44.h
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -107,6 +107,7 @@ Revision 1.12 2004/03/25 14:57:49 ponchio
#include <vcg/space/point3.h> #include <vcg/space/point3.h>
#include <vcg/space/point4.h> #include <vcg/space/point4.h>
#include <vector> #include <vector>
#include <iostream>
namespace vcg { namespace vcg {
@ -120,13 +121,13 @@ Opengl stores matrix in column-major order. That is, the matrix is stored as:
a2 a6 a10 a14 a2 a6 a10 a14
a3 a7 a11 a15 a3 a7 a11 a15
Usually in opengl (see opengl specs) vectors are 'column' vectors Usually in opengl (see opengl specs) vectors are 'column' vectors
so usually matrix are PRE-multiplied for a vector. so usually matrix are PRE-multiplied for a vector.
So the command glTranslate generate a matrix that So the command glTranslate generate a matrix that
is ready to be premultipled for a vector: is ready to be premultipled for a vector:
1 0 0 tx 1 0 0 tx
0 1 0 ty 0 1 0 ty
0 0 1 tz 0 0 1 tz
0 0 0 1 0 0 0 1
@ -138,7 +139,7 @@ Matrix44 stores matrix in row-major order i.e.
a12 a13 a14 a15 a12 a13 a14 a15
So for the use of that matrix in opengl with their supposed meaning you have to transpose them before feeding to glMultMatrix. So for the use of that matrix in opengl with their supposed meaning you have to transpose them before feeding to glMultMatrix.
This mechanism is hidden by the templated function defined in wrap/gl/math.h; This mechanism is hidden by the templated function defined in wrap/gl/math.h;
If your machine has the ARB_transpose_matrix extension it will use the appropriate; If your machine has the ARB_transpose_matrix extension it will use the appropriate;
The various gl-like command SetRotate, SetTranslate assume that you are making matrix The various gl-like command SetRotate, SetTranslate assume that you are making matrix
for 'column' vectors. for 'column' vectors.
@ -150,7 +151,7 @@ class Matrix44Diag:public Point4<S>{
public: public:
/** @name Matrix33 /** @name Matrix33
Class Matrix33Diag. Class Matrix33Diag.
This is the class for definition of a diagonal matrix 4x4. This is the class for definition of a diagonal matrix 4x4.
@param S (Templete Parameter) Specifies the ScalarType field. @param S (Templete Parameter) Specifies the ScalarType field.
*/ */
Matrix44Diag(const S & p0,const S & p1,const S & p2,const S & p3):Point4<S>(p0,p1,p2,p3){}; Matrix44Diag(const S & p0,const S & p1,const S & p2,const S & p3):Point4<S>(p0,p1,p2,p3){};
@ -160,11 +161,11 @@ public:
/** This class represent a 4x4 matrix. T is the kind of element in the matrix. /** This class represent a 4x4 matrix. T is the kind of element in the matrix.
*/ */
template <class T> class Matrix44 { template <class T> class Matrix44 {
protected: protected:
T _a[16]; T _a[16];
public: public:
typedef T ScalarType; typedef T ScalarType;
///@{ ///@{
@ -172,7 +173,7 @@ public:
/** $name Constructors /** $name Constructors
* No automatic casting and default constructor is empty * No automatic casting and default constructor is empty
*/ */
Matrix44() {}; Matrix44() {};
~Matrix44() {}; ~Matrix44() {};
Matrix44(const Matrix44 &m); Matrix44(const Matrix44 &m);
Matrix44(const T v[]); Matrix44(const T v[]);
@ -195,7 +196,7 @@ public:
//const T &operator[](const int i) const; //const T &operator[](const int i) const;
T *V(); T *V();
const T *V() const ; const T *V() const ;
T *operator[](const int i); T *operator[](const int i);
const T *operator[](const int i) const; const T *operator[](const int i) const;
@ -216,7 +217,7 @@ public:
return Point4<T>(ElementAt(i,0),ElementAt(i,1),ElementAt(i,2),ElementAt(i,3)); return Point4<T>(ElementAt(i,0),ElementAt(i,1),ElementAt(i,2),ElementAt(i,3));
// return *((Point4<T>*)(&_a[i<<2])); alternativa forse + efficiente // return *((Point4<T>*)(&_a[i<<2])); alternativa forse + efficiente
} }
Point3<T> GetRow3(const int& i)const{ Point3<T> GetRow3(const int& i)const{
assert(i>=0 && i<4); assert(i>=0 && i<4);
return Point3<T>(ElementAt(i,0),ElementAt(i,1),ElementAt(i,2)); return Point3<T>(ElementAt(i,0),ElementAt(i,1),ElementAt(i,2));
@ -227,18 +228,18 @@ public:
Matrix44 operator-(const Matrix44 &m) const; Matrix44 operator-(const Matrix44 &m) const;
Matrix44 operator*(const Matrix44 &m) const; Matrix44 operator*(const Matrix44 &m) const;
Matrix44 operator*(const Matrix44Diag<T> &m) const; Matrix44 operator*(const Matrix44Diag<T> &m) const;
Point4<T> operator*(const Point4<T> &v) const; Point4<T> operator*(const Point4<T> &v) const;
bool operator==(const Matrix44 &m) const; bool operator==(const Matrix44 &m) const;
bool operator!= (const Matrix44 &m) const; bool operator!= (const Matrix44 &m) const;
Matrix44 operator-() const; Matrix44 operator-() const;
Matrix44 operator*(const T k) const; Matrix44 operator*(const T k) const;
void operator+=(const Matrix44 &m); void operator+=(const Matrix44 &m);
void operator-=(const Matrix44 &m); void operator-=(const Matrix44 &m);
void operator*=( const Matrix44 & m ); void operator*=( const Matrix44 & m );
void operator*=( const T k ); void operator*=( const T k );
template <class Matrix44Type> template <class Matrix44Type>
void ToMatrix(Matrix44Type & m) const {for(int i = 0; i < 16; i++) m.V()[i]=V()[i];} void ToMatrix(Matrix44Type & m) const {for(int i = 0; i < 16; i++) m.V()[i]=V()[i];}
@ -250,31 +251,31 @@ public:
void SetZero(); void SetZero();
void SetIdentity(); void SetIdentity();
void SetDiagonal(const T k); void SetDiagonal(const T k);
Matrix44 &SetScale(const T sx, const T sy, const T sz); Matrix44 &SetScale(const T sx, const T sy, const T sz);
Matrix44 &SetScale(const Point3<T> &t); Matrix44 &SetScale(const Point3<T> &t);
Matrix44 &SetTranslate(const Point3<T> &t); Matrix44 &SetTranslate(const Point3<T> &t);
Matrix44 &SetTranslate(const T sx, const T sy, const T sz); Matrix44 &SetTranslate(const T sx, const T sy, const T sz);
Matrix44 &SetShearXY(const T sz); Matrix44 &SetShearXY(const T sz);
Matrix44 &SetShearXZ(const T sy); Matrix44 &SetShearXZ(const T sy);
Matrix44 &SetShearYZ(const T sx); Matrix44 &SetShearYZ(const T sx);
///use radiants for angle. ///use radiants for angle.
Matrix44 &SetRotateDeg(T AngleDeg, const Point3<T> & axis); Matrix44 &SetRotateDeg(T AngleDeg, const Point3<T> & axis);
Matrix44 &SetRotateRad(T AngleRad, const Point3<T> & axis); Matrix44 &SetRotateRad(T AngleRad, const Point3<T> & axis);
T Determinant() const; T Determinant() const;
template <class Q> void Import(const Matrix44<Q> &m) { template <class Q> void Import(const Matrix44<Q> &m) {
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
_a[i] = (T)(m.V()[i]); _a[i] = (T)(m.V()[i]);
} }
template <class Q> template <class Q>
static inline Matrix44 Construct( const Matrix44<Q> & b ) static inline Matrix44 Construct( const Matrix44<Q> & b )
{ {
Matrix44<T> tmp; tmp.FromMatrix(b); Matrix44<T> tmp; tmp.FromMatrix(b);
return tmp; return tmp;
} }
static inline const Matrix44 &Identity( ) static inline const Matrix44 &Identity( )
{ {
static Matrix44<T> tmp; tmp.SetIdentity(); static Matrix44<T> tmp; tmp.SetIdentity();
@ -291,6 +292,18 @@ public:
// for the transistion to eigen // for the transistion to eigen
const Matrix44& eval() { return *this; } const Matrix44& eval() { return *this; }
void print() {
unsigned int i, j, p;
for (i=0, p=0; i<4; i++, p+=4)
{
std::cout << "[\t";
for (j=0; j<4; j++)
std::cout << _a[p+j] << "\t";
std::cout << "]\n";
}
std::cout << "\n";
}
}; };
@ -299,7 +312,7 @@ public:
template <class T> class LinearSolve: public Matrix44<T> { template <class T> class LinearSolve: public Matrix44<T> {
public: public:
LinearSolve(const Matrix44<T> &m); LinearSolve(const Matrix44<T> &m);
Point4<T> Solve(const Point4<T> &b); // solve A <20> x = b Point4<T> Solve(const Point4<T> &b); // solve A <20> x = b
///If you need to solve some equation you can use this function instead of Matrix44 one for speed. ///If you need to solve some equation you can use this function instead of Matrix44 one for speed.
T Determinant() const; T Determinant() const;
protected: protected:
@ -313,13 +326,13 @@ protected:
/*** Postmultiply */ /*** Postmultiply */
//template <class T> Point3<T> operator*(const Point3<T> &p, const Matrix44<T> &m); //template <class T> Point3<T> operator*(const Point3<T> &p, const Matrix44<T> &m);
///Premultiply ///Premultiply
template <class T> Point3<T> operator*(const Matrix44<T> &m, const Point3<T> &p); template <class T> Point3<T> operator*(const Matrix44<T> &m, const Point3<T> &p);
template <class T> Matrix44<T> &Transpose(Matrix44<T> &m); template <class T> Matrix44<T> &Transpose(Matrix44<T> &m);
//return NULL matrix if not invertible //return NULL matrix if not invertible
template <class T> Matrix44<T> &Invert(Matrix44<T> &m); template <class T> Matrix44<T> &Invert(Matrix44<T> &m);
template <class T> Matrix44<T> Inverse(const Matrix44<T> &m); template <class T> Matrix44<T> Inverse(const Matrix44<T> &m);
typedef Matrix44<short> Matrix44s; typedef Matrix44<short> Matrix44s;
typedef Matrix44<int> Matrix44i; typedef Matrix44<int> Matrix44i;
@ -329,11 +342,11 @@ typedef Matrix44<double> Matrix44d;
template <class T> Matrix44<T>::Matrix44(const Matrix44<T> &m) { template <class T> Matrix44<T>::Matrix44(const Matrix44<T> &m) {
memcpy((T *)_a, (T *)m._a, 16 * sizeof(T)); memcpy((T *)_a, (T *)m._a, 16 * sizeof(T));
} }
template <class T> Matrix44<T>::Matrix44(const T v[]) { template <class T> Matrix44<T>::Matrix44(const T v[]) {
memcpy((T *)_a, v, 16 * sizeof(T)); memcpy((T *)_a, v, 16 * sizeof(T));
} }
template <class T> T &Matrix44<T>::ElementAt(const int row, const int col) { template <class T> T &Matrix44<T>::ElementAt(const int row, const int col) {
@ -356,7 +369,7 @@ template <class T> T Matrix44<T>::ElementAt(const int row, const int col) const
//template <class T> const T &Matrix44<T>::operator[](const int i) const { //template <class T> const T &Matrix44<T>::operator[](const int i) const {
// assert(i >= 0 && i < 16); // assert(i >= 0 && i < 16);
// return ((T *)_a)[i]; // return ((T *)_a)[i];
//} //}
template <class T> T *Matrix44<T>::operator[](const int i) { template <class T> T *Matrix44<T>::operator[](const int i) {
assert(i >= 0 && i < 4); assert(i >= 0 && i < 4);
return _a+i*4; return _a+i*4;
@ -366,26 +379,26 @@ template <class T> const T *Matrix44<T>::operator[](const int i) const {
assert(i >= 0 && i < 4); assert(i >= 0 && i < 4);
return _a+i*4; return _a+i*4;
} }
template <class T> T *Matrix44<T>::V() { return _a;} template <class T> T *Matrix44<T>::V() { return _a;}
template <class T> const T *Matrix44<T>::V() const { return _a;} template <class T> const T *Matrix44<T>::V() const { return _a;}
template <class T> Matrix44<T> Matrix44<T>::operator+(const Matrix44 &m) const { template <class T> Matrix44<T> Matrix44<T>::operator+(const Matrix44 &m) const {
Matrix44<T> ret; Matrix44<T> ret;
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
ret.V()[i] = V()[i] + m.V()[i]; ret.V()[i] = V()[i] + m.V()[i];
return ret; return ret;
} }
template <class T> Matrix44<T> Matrix44<T>::operator-(const Matrix44 &m) const { template <class T> Matrix44<T> Matrix44<T>::operator-(const Matrix44 &m) const {
Matrix44<T> ret; Matrix44<T> ret;
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
ret.V()[i] = V()[i] - m.V()[i]; ret.V()[i] = V()[i] - m.V()[i];
return ret; return ret;
} }
template <class T> Matrix44<T> Matrix44<T>::operator*(const Matrix44 &m) const { template <class T> Matrix44<T> Matrix44<T>::operator*(const Matrix44 &m) const {
Matrix44 ret; Matrix44 ret;
for(int i = 0; i < 4; i++) for(int i = 0; i < 4; i++)
for(int j = 0; j < 4; j++) { for(int j = 0; j < 4; j++) {
T t = 0.0; T t = 0.0;
@ -397,15 +410,15 @@ template <class T> Matrix44<T> Matrix44<T>::operator*(const Matrix44 &m) const {
} }
template <class T> Matrix44<T> Matrix44<T>::operator*(const Matrix44Diag<T> &m) const { template <class T> Matrix44<T> Matrix44<T>::operator*(const Matrix44Diag<T> &m) const {
Matrix44 ret = (*this); Matrix44 ret = (*this);
for(int i = 0; i < 4; ++i) for(int i = 0; i < 4; ++i)
for(int j = 0; j < 4; ++j) for(int j = 0; j < 4; ++j)
ret[i][j]*=m[i]; ret[i][j]*=m[i];
return ret; return ret;
} }
template <class T> Point4<T> Matrix44<T>::operator*(const Point4<T> &v) const { template <class T> Point4<T> Matrix44<T>::operator*(const Point4<T> &v) const {
Point4<T> ret; Point4<T> ret;
for(int i = 0; i < 4; i++){ for(int i = 0; i < 4; i++){
T t = 0.0; T t = 0.0;
for(int k = 0; k < 4; k++) for(int k = 0; k < 4; k++)
@ -418,14 +431,14 @@ template <class T> Point4<T> Matrix44<T>::operator*(const Point4<T> &v) const {
template <class T> bool Matrix44<T>::operator==(const Matrix44 &m) const { template <class T> bool Matrix44<T>::operator==(const Matrix44 &m) const {
for(int i = 0; i < 4; ++i) for(int i = 0; i < 4; ++i)
for(int j = 0; j < 4; ++j) for(int j = 0; j < 4; ++j)
if(ElementAt(i,j) != m.ElementAt(i,j)) if(ElementAt(i,j) != m.ElementAt(i,j))
return false; return false;
return true; return true;
} }
template <class T> bool Matrix44<T>::operator!=(const Matrix44 &m) const { template <class T> bool Matrix44<T>::operator!=(const Matrix44 &m) const {
for(int i = 0; i < 4; ++i) for(int i = 0; i < 4; ++i)
for(int j = 0; j < 4; ++j) for(int j = 0; j < 4; ++j)
if(ElementAt(i,j) != m.ElementAt(i,j)) if(ElementAt(i,j) != m.ElementAt(i,j))
return true; return true;
return false; return false;
@ -447,7 +460,7 @@ template <class T> Matrix44<T> Matrix44<T>::operator*(const T k) const {
template <class T> void Matrix44<T>::operator+=(const Matrix44 &m) { template <class T> void Matrix44<T>::operator+=(const Matrix44 &m) {
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
V()[i] += m.V()[i]; V()[i] += m.V()[i];
} }
template <class T> void Matrix44<T>::operator-=(const Matrix44 &m) { template <class T> void Matrix44<T>::operator-=(const Matrix44 &m) {
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
@ -455,7 +468,7 @@ template <class T> void Matrix44<T>::operator-=(const Matrix44 &m) {
} }
template <class T> void Matrix44<T>::operator*=( const Matrix44 & m ) { template <class T> void Matrix44<T>::operator*=( const Matrix44 & m ) {
*this = *this *m; *this = *this *m;
/*for(int i = 0; i < 4; i++) { //sbagliato /*for(int i = 0; i < 4; i++) { //sbagliato
Point4<T> t(0, 0, 0, 0); Point4<T> t(0, 0, 0, 0);
for(int k = 0; k < 4; k++) { for(int k = 0; k < 4; k++) {
@ -487,7 +500,7 @@ void Matrix44<T>::ToEulerAngles(T &alpha, T &beta, T &gamma)
gamma = atan2(ElementAt(0,1), ElementAt(1,1)); gamma = atan2(ElementAt(0,1), ElementAt(1,1));
} }
template <class T> template <class T>
void Matrix44<T>::FromEulerAngles(T alpha, T beta, T gamma) void Matrix44<T>::FromEulerAngles(T alpha, T beta, T gamma)
{ {
this->SetZero(); this->SetZero();
@ -499,18 +512,18 @@ void Matrix44<T>::FromEulerAngles(T alpha, T beta, T gamma)
T sinbeta = sin(beta); T sinbeta = sin(beta);
T singamma = sin(gamma); T singamma = sin(gamma);
ElementAt(0,0) = cosbeta * cosgamma; ElementAt(0,0) = cosbeta * cosgamma;
ElementAt(1,0) = -cosalpha * singamma + sinalpha * sinbeta * cosgamma; ElementAt(1,0) = -cosalpha * singamma + sinalpha * sinbeta * cosgamma;
ElementAt(2,0) = sinalpha * singamma + cosalpha * sinbeta * cosgamma; ElementAt(2,0) = sinalpha * singamma + cosalpha * sinbeta * cosgamma;
ElementAt(0,1) = cosbeta * singamma; ElementAt(0,1) = cosbeta * singamma;
ElementAt(1,1) = cosalpha * cosgamma + sinalpha * sinbeta * singamma; ElementAt(1,1) = cosalpha * cosgamma + sinalpha * sinbeta * singamma;
ElementAt(2,1) = -sinalpha * cosgamma + cosalpha * sinbeta * singamma; ElementAt(2,1) = -sinalpha * cosgamma + cosalpha * sinbeta * singamma;
ElementAt(0,2) = -sinbeta; ElementAt(0,2) = -sinbeta;
ElementAt(1,2) = sinalpha * cosbeta; ElementAt(1,2) = sinalpha * cosbeta;
ElementAt(2,2) = cosalpha * cosbeta; ElementAt(2,2) = cosalpha * cosbeta;
ElementAt(3,3) = 1; ElementAt(3,3) = 1;
} }
@ -518,8 +531,8 @@ template <class T> void Matrix44<T>::SetZero() {
memset((T *)_a, 0, 16 * sizeof(T)); memset((T *)_a, 0, 16 * sizeof(T));
} }
template <class T> void Matrix44<T>::SetIdentity() { template <class T> void Matrix44<T>::SetIdentity() {
SetDiagonal(1); SetDiagonal(1);
} }
template <class T> void Matrix44<T>::SetDiagonal(const T k) { template <class T> void Matrix44<T>::SetDiagonal(const T k) {
@ -527,7 +540,7 @@ template <class T> void Matrix44<T>::SetDiagonal(const T k) {
ElementAt(0, 0) = k; ElementAt(0, 0) = k;
ElementAt(1, 1) = k; ElementAt(1, 1) = k;
ElementAt(2, 2) = k; ElementAt(2, 2) = k;
ElementAt(3, 3) = 1; ElementAt(3, 3) = 1;
} }
template <class T> Matrix44<T> &Matrix44<T>::SetScale(const Point3<T> &t) { template <class T> Matrix44<T> &Matrix44<T>::SetScale(const Point3<T> &t) {
@ -555,15 +568,15 @@ template <class T> Matrix44<T> &Matrix44<T>::SetTranslate(const T tx, const T ty
return *this; return *this;
} }
template <class T> Matrix44<T> &Matrix44<T>::SetRotateDeg(T AngleDeg, const Point3<T> & axis) { template <class T> Matrix44<T> &Matrix44<T>::SetRotateDeg(T AngleDeg, const Point3<T> & axis) {
return SetRotateRad(math::ToRad(AngleDeg),axis); return SetRotateRad(math::ToRad(AngleDeg),axis);
} }
template <class T> Matrix44<T> &Matrix44<T>::SetRotateRad(T AngleRad, const Point3<T> & axis) { template <class T> Matrix44<T> &Matrix44<T>::SetRotateRad(T AngleRad, const Point3<T> & axis) {
//angle = angle*(T)3.14159265358979323846/180; e' in radianti! //angle = angle*(T)3.14159265358979323846/180; e' in radianti!
T c = math::Cos(AngleRad); T c = math::Cos(AngleRad);
T s = math::Sin(AngleRad); T s = math::Sin(AngleRad);
T q = 1-c; T q = 1-c;
Point3<T> t = axis; Point3<T> t = axis;
t.Normalize(); t.Normalize();
ElementAt(0,0) = t[0]*t[0]*q + c; ElementAt(0,0) = t[0]*t[0]*q + c;
@ -579,14 +592,14 @@ template <class T> Matrix44<T> &Matrix44<T>::SetRotateRad(T AngleRad, const Poin
ElementAt(2,2) = t[2]*t[2]*q +c; ElementAt(2,2) = t[2]*t[2]*q +c;
ElementAt(2,3) = 0; ElementAt(2,3) = 0;
ElementAt(3,0) = 0; ElementAt(3,0) = 0;
ElementAt(3,1) = 0; ElementAt(3,1) = 0;
ElementAt(3,2) = 0; ElementAt(3,2) = 0;
ElementAt(3,3) = 1; ElementAt(3,3) = 1;
return *this; return *this;
} }
/* Shear Matrixes /* Shear Matrixes
XY XY
1 k 0 0 x x+ky 1 k 0 0 x x+ky
0 1 0 0 y y 0 1 0 0 y y
0 0 1 0 z z 0 0 1 0 z z
@ -604,19 +617,19 @@ template <class T> Matrix44<T> &Matrix44<T>::SetRotateRad(T AngleRad, const Poin
*/ */
template <class T> Matrix44<T> & Matrix44<T>:: SetShearXY( const T sh) {// shear the X coordinate as the Y coordinate change template <class T> Matrix44<T> & Matrix44<T>:: SetShearXY( const T sh) {// shear the X coordinate as the Y coordinate change
SetIdentity(); SetIdentity();
ElementAt(0,1) = sh; ElementAt(0,1) = sh;
return *this; return *this;
} }
template <class T> Matrix44<T> & Matrix44<T>:: SetShearXZ( const T sh) {// shear the X coordinate as the Z coordinate change template <class T> Matrix44<T> & Matrix44<T>:: SetShearXZ( const T sh) {// shear the X coordinate as the Z coordinate change
SetIdentity(); SetIdentity();
ElementAt(0,2) = sh; ElementAt(0,2) = sh;
return *this; return *this;
} }
template <class T> Matrix44<T> &Matrix44<T>:: SetShearYZ( const T sh) {// shear the Y coordinate as the Z coordinate change template <class T> Matrix44<T> &Matrix44<T>:: SetShearYZ( const T sh) {// shear the Y coordinate as the Z coordinate change
SetIdentity(); SetIdentity();
ElementAt(1,2) = sh; ElementAt(1,2) = sh;
return *this; return *this;
@ -625,7 +638,7 @@ template <class T> Matrix44<T> &Matrix44<T>::SetRotateRad(T AngleRad, const Poin
/* /*
Given a non singular, non projective matrix (e.g. with the last row equal to [0,0,0,1] ) Given a non singular, non projective matrix (e.g. with the last row equal to [0,0,0,1] )
This procedure decompose it in a sequence of This procedure decompose it in a sequence of
Scale,Shear,Rotation e Translation Scale,Shear,Rotation e Translation
- ScaleV and Tranv are obiviously scaling and translation. - ScaleV and Tranv are obiviously scaling and translation.
@ -635,7 +648,7 @@ This procedure decompose it in a sequence of
The input matrix is modified leaving inside it a simple roto translation. The input matrix is modified leaving inside it a simple roto translation.
To obtain the original matrix the above transformation have to be applied in the strict following way: To obtain the original matrix the above transformation have to be applied in the strict following way:
OriginalMatrix = Trn * Rtx*Rty*Rtz * ShearYZ*ShearXZ*ShearXY * Scl OriginalMatrix = Trn * Rtx*Rty*Rtz * ShearYZ*ShearXZ*ShearXY * Scl
Example Code: Example Code:
@ -655,7 +668,7 @@ double srv() { return (double(rand()%40)-20)/2.0; } // small random value
Matrix44d Rty; Rty.SetRotate(math::ToRad(RtV[1]),Point3d(0,1,0)); Matrix44d Rty; Rty.SetRotate(math::ToRad(RtV[1]),Point3d(0,1,0));
Matrix44d Rtz; Rtz.SetRotate(math::ToRad(RtV[2]),Point3d(0,0,1)); Matrix44d Rtz; Rtz.SetRotate(math::ToRad(RtV[2]),Point3d(0,0,1));
Matrix44d Trn; Trn.SetTranslate(TrV); Matrix44d Trn; Trn.SetTranslate(TrV);
Matrix44d StartM = Trn * Rtx*Rty*Rtz * Syz*Sxz*Sxy *Scl; Matrix44d StartM = Trn * Rtx*Rty*Rtz * Syz*Sxz*Sxy *Scl;
Matrix44d ResultM=StartM; Matrix44d ResultM=StartM;
Decompose(ResultM,ScVOut,ShVOut,RtVOut,TrVOut); Decompose(ResultM,ScVOut,ShVOut,RtVOut,TrVOut);
@ -670,75 +683,75 @@ double srv() { return (double(rand()%40)-20)/2.0; } // small random value
Trn.SetTranslate(TrVOut); Trn.SetTranslate(TrVOut);
// Now Rebuild is equal to StartM // Now Rebuild is equal to StartM
Matrix44d RebuildM = Trn * Rtx*Rty*Rtz * Syz*Sxz*Sxy * Scl ; Matrix44d RebuildM = Trn * Rtx*Rty*Rtz * Syz*Sxz*Sxy * Scl ;
*/ */
template <class T> template <class T>
bool Decompose(Matrix44<T> &M, Point3<T> &ScaleV, Point3<T> &ShearV, Point3<T> &RotV,Point3<T> &TranV) bool Decompose(Matrix44<T> &M, Point3<T> &ScaleV, Point3<T> &ShearV, Point3<T> &RotV,Point3<T> &TranV)
{ {
if(!(M[3][0]==0 && M[3][1]==0 && M[3][2]==0 && M[3][3]==1) ) // the matrix is projective if(!(M[3][0]==0 && M[3][1]==0 && M[3][2]==0 && M[3][3]==1) ) // the matrix is projective
return false; return false;
if(math::Abs(M.Determinant())<1e-10) return false; // matrix should be at least invertible... if(math::Abs(M.Determinant())<1e-10) return false; // matrix should be at least invertible...
// First Step recover the traslation // First Step recover the traslation
TranV=M.GetColumn3(3); TranV=M.GetColumn3(3);
// Second Step Recover Scale and Shearing interleaved // Second Step Recover Scale and Shearing interleaved
ScaleV[0]=Norm(M.GetColumn3(0)); ScaleV[0]=Norm(M.GetColumn3(0));
Point3<T> R[3]; Point3<T> R[3];
R[0]=M.GetColumn3(0); R[0]=M.GetColumn3(0);
R[0].Normalize(); R[0].Normalize();
ShearV[0]=R[0]*M.GetColumn3(1); // xy shearing ShearV[0]=R[0]*M.GetColumn3(1); // xy shearing
R[1]= M.GetColumn3(1)-R[0]*ShearV[0]; R[1]= M.GetColumn3(1)-R[0]*ShearV[0];
assert(math::Abs(R[1]*R[0])<1e-10); assert(math::Abs(R[1]*R[0])<1e-10);
ScaleV[1]=Norm(R[1]); // y scaling ScaleV[1]=Norm(R[1]); // y scaling
R[1]=R[1]/ScaleV[1]; R[1]=R[1]/ScaleV[1];
ShearV[0]=ShearV[0]/ScaleV[1]; ShearV[0]=ShearV[0]/ScaleV[1];
ShearV[1]=R[0]*M.GetColumn3(2); // xz shearing ShearV[1]=R[0]*M.GetColumn3(2); // xz shearing
R[2]= M.GetColumn3(2)-R[0]*ShearV[1]; R[2]= M.GetColumn3(2)-R[0]*ShearV[1];
assert(math::Abs(R[2]*R[0])<1e-10); assert(math::Abs(R[2]*R[0])<1e-10);
R[2] = R[2]-R[1]*(R[2]*R[1]); R[2] = R[2]-R[1]*(R[2]*R[1]);
assert(math::Abs(R[2]*R[1])<1e-10); assert(math::Abs(R[2]*R[1])<1e-10);
assert(math::Abs(R[2]*R[0])<1e-10); assert(math::Abs(R[2]*R[0])<1e-10);
ScaleV[2]=Norm(R[2]); ScaleV[2]=Norm(R[2]);
ShearV[1]=ShearV[1]/ScaleV[2]; ShearV[1]=ShearV[1]/ScaleV[2];
R[2]=R[2]/ScaleV[2]; R[2]=R[2]/ScaleV[2];
assert(math::Abs(R[2]*R[1])<1e-10); assert(math::Abs(R[2]*R[1])<1e-10);
assert(math::Abs(R[2]*R[0])<1e-10); assert(math::Abs(R[2]*R[0])<1e-10);
ShearV[2]=R[1]*M.GetColumn3(2); // yz shearing ShearV[2]=R[1]*M.GetColumn3(2); // yz shearing
ShearV[2]=ShearV[2]/ScaleV[2]; ShearV[2]=ShearV[2]/ScaleV[2];
int i,j; int i,j;
for(i=0;i<3;++i) for(i=0;i<3;++i)
for(j=0;j<3;++j) for(j=0;j<3;++j)
M[i][j]=R[j][i]; M[i][j]=R[j][i];
// Third and last step: Recover the rotation // Third and last step: Recover the rotation
//now the matrix should be a pure rotation matrix so its determinant is +-1 //now the matrix should be a pure rotation matrix so its determinant is +-1
double det=M.Determinant(); double det=M.Determinant();
if(math::Abs(det)<1e-10) return false; // matrix should be at least invertible... if(math::Abs(det)<1e-10) return false; // matrix should be at least invertible...
assert(math::Abs(math::Abs(det)-1.0)<1e-10); // it should be +-1... assert(math::Abs(math::Abs(det)-1.0)<1e-10); // it should be +-1...
if(det<0) { if(det<0) {
ScaleV *= -1; ScaleV *= -1;
M *= -1; M *= -1;
} }
double alpha,beta,gamma; // rotations around the x,y and z axis double alpha,beta,gamma; // rotations around the x,y and z axis
beta=asin( M[0][2]); beta=asin( M[0][2]);
double cosbeta=cos(beta); double cosbeta=cos(beta);
if(math::Abs(cosbeta) > 1e-5) if(math::Abs(cosbeta) > 1e-5)
{ {
alpha=asin(-M[1][2]/cosbeta); alpha=asin(-M[1][2]/cosbeta);
if((M[2][2]/cosbeta) < 0 ) alpha=M_PI-alpha; if((M[2][2]/cosbeta) < 0 ) alpha=M_PI-alpha;
gamma=asin(-M[0][1]/cosbeta); gamma=asin(-M[0][1]/cosbeta);
if((M[0][0]/cosbeta)<0) gamma = M_PI-gamma; if((M[0][0]/cosbeta)<0) gamma = M_PI-gamma;
} }
else else
{ {
alpha=asin(-M[1][0]); alpha=asin(-M[1][0]);
if(M[1][1]<0) alpha=M_PI-alpha; if(M[1][1]<0) alpha=M_PI-alpha;
gamma=0; gamma=0;
@ -754,7 +767,7 @@ bool Decompose(Matrix44<T> &M, Point3<T> &ScaleV, Point3<T> &ShearV, Point3<T> &
template <class T> T Matrix44<T>::Determinant() const { template <class T> T Matrix44<T>::Determinant() const {
LinearSolve<T> solve(*this); LinearSolve<T> solve(*this);
return solve.Determinant(); return solve.Determinant();
} }
@ -785,45 +798,45 @@ template <class T> Point3<T> operator*(const Matrix44<T> &m, const Point3<T> &p)
template <class T> Matrix44<T> &Transpose(Matrix44<T> &m) { template <class T> Matrix44<T> &Transpose(Matrix44<T> &m) {
for(int i = 1; i < 4; i++) for(int i = 1; i < 4; i++)
for(int j = 0; j < i; j++) { for(int j = 0; j < i; j++) {
math::Swap(m.ElementAt(i, j), m.ElementAt(j, i)); math::Swap(m.ElementAt(i, j), m.ElementAt(j, i));
} }
return m; return m;
} }
/* /*
To invert a matrix you can To invert a matrix you can
either invert the matrix inplace calling either invert the matrix inplace calling
vcg::Invert(yourMatrix); vcg::Invert(yourMatrix);
or get the inverse matrix of a given matrix without touching it: or get the inverse matrix of a given matrix without touching it:
invertedMatrix = vcg::Inverse(untouchedMatrix); invertedMatrix = vcg::Inverse(untouchedMatrix);
*/ */
template <class T> Matrix44<T> & Invert(Matrix44<T> &m) { template <class T> Matrix44<T> & Invert(Matrix44<T> &m) {
LinearSolve<T> solve(m); LinearSolve<T> solve(m);
for(int j = 0; j < 4; j++) { //Find inverse by columns. for(int j = 0; j < 4; j++) { //Find inverse by columns.
Point4<T> col(0, 0, 0, 0); Point4<T> col(0, 0, 0, 0);
col[j] = 1.0; col[j] = 1.0;
col = solve.Solve(col); col = solve.Solve(col);
for(int i = 0; i < 4; i++) for(int i = 0; i < 4; i++)
m.ElementAt(i, j) = col[i]; m.ElementAt(i, j) = col[i];
} }
return m; return m;
} }
template <class T> Matrix44<T> Inverse(const Matrix44<T> &m) { template <class T> Matrix44<T> Inverse(const Matrix44<T> &m) {
LinearSolve<T> solve(m); LinearSolve<T> solve(m);
Matrix44<T> res; Matrix44<T> res;
for(int j = 0; j < 4; j++) { //Find inverse by columns. for(int j = 0; j < 4; j++) { //Find inverse by columns.
Point4<T> col(0, 0, 0, 0); Point4<T> col(0, 0, 0, 0);
col[j] = 1.0; col[j] = 1.0;
col = solve.Solve(col); col = solve.Solve(col);
for(int i = 0; i < 4; i++) for(int i = 0; i < 4; i++)
res.ElementAt(i, j) = col[i]; res.ElementAt(i, j) = col[i];
} }
return res; return res;
} }
@ -842,8 +855,8 @@ template <class T> LinearSolve<T>::LinearSolve(const Matrix44<T> &m): Matrix44<T
template <class T> T LinearSolve<T>::Determinant() const { template <class T> T LinearSolve<T>::Determinant() const {
T det = d; T det = d;
for(int j = 0; j < 4; j++) for(int j = 0; j < 4; j++)
det *= this-> ElementAt(j, j); det *= this-> ElementAt(j, j);
return det; return det;
} }
@ -852,18 +865,18 @@ template <class T> T LinearSolve<T>::Determinant() const {
d is +or -1 depeneing of row permutation even or odd.*/ d is +or -1 depeneing of row permutation even or odd.*/
#define TINY 1e-100 #define TINY 1e-100
template <class T> bool LinearSolve<T>::Decompose() { template <class T> bool LinearSolve<T>::Decompose() {
/* Matrix44<T> A; /* Matrix44<T> A;
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
A[i] = operator[](i); A[i] = operator[](i);
SetIdentity(); SetIdentity();
Point4<T> scale; Point4<T> scale;
// Set scale factor, scale(i) = max( |a(i,j)| ), for each row // Set scale factor, scale(i) = max( |a(i,j)| ), for each row
for(int i = 0; i < 4; i++ ) { for(int i = 0; i < 4; i++ ) {
index[i] = i; // Initialize row index list index[i] = i; // Initialize row index list
T scalemax = (T)0.0; T scalemax = (T)0.0;
for(int j = 0; j < 4; j++) for(int j = 0; j < 4; j++)
scalemax = (scalemax > math::Abs(A.ElementAt(i,j))) ? scalemax : math::Abs(A.ElementAt(i,j)); scalemax = (scalemax > math::Abs(A.ElementAt(i,j))) ? scalemax : math::Abs(A.ElementAt(i,j));
scale[i] = scalemax; scale[i] = scalemax;
} }
@ -895,23 +908,23 @@ template <class T> bool LinearSolve<T>::Decompose() {
for(int j=k+1; j < 4; j++ ) for(int j=k+1; j < 4; j++ )
A.ElementAt(index[i],j) -= coeff*A.ElementAt(indexJ,j); A.ElementAt(index[i],j) -= coeff*A.ElementAt(indexJ,j);
A.ElementAt(index[i],k) = coeff; A.ElementAt(index[i],k) = coeff;
//for( j=1; j< 4; j++ ) //for( j=1; j< 4; j++ )
// ElementAt(index[i],j) -= A.ElementAt(index[i], k)*ElementAt(indexJ, j); // ElementAt(index[i],j) -= A.ElementAt(index[i], k)*ElementAt(indexJ, j);
} }
} }
for(int i = 0; i < 16; i++) for(int i = 0; i < 16; i++)
operator[](i) = A[i]; operator[](i) = A[i];
d = signDet; d = signDet;
// this = A; // this = A;
return true; */ return true; */
d = 1; //no permutation still d = 1; //no permutation still
T scaling[4]; T scaling[4];
int i,j,k; int i,j,k;
//Saving the scvaling information per row //Saving the scvaling information per row
for(i = 0; i < 4; i++) { for(i = 0; i < 4; i++) {
T largest = 0.0; T largest = 0.0;
for(j = 0; j < 4; j++) { for(j = 0; j < 4; j++) {
T t = math::Abs(this->ElementAt(i, j)); T t = math::Abs(this->ElementAt(i, j));
@ -920,70 +933,70 @@ template <class T> bool LinearSolve<T>::Decompose() {
if (largest == 0.0) { //oooppps there is a zero row! if (largest == 0.0) { //oooppps there is a zero row!
return false; return false;
} }
scaling[i] = (T)1.0 / largest; scaling[i] = (T)1.0 / largest;
} }
int imax = 0; int imax = 0;
for(j = 0; j < 4; j++) { for(j = 0; j < 4; j++) {
for(i = 0; i < j; i++) { for(i = 0; i < j; i++) {
T sum = this->ElementAt(i,j); T sum = this->ElementAt(i,j);
for(int k = 0; k < i; k++) for(int k = 0; k < i; k++)
sum -= this->ElementAt(i,k)*this->ElementAt(k,j); sum -= this->ElementAt(i,k)*this->ElementAt(k,j);
this->ElementAt(i,j) = sum; this->ElementAt(i,j) = sum;
} }
T largest = 0.0; T largest = 0.0;
for(i = j; i < 4; i++) { for(i = j; i < 4; i++) {
T sum = this->ElementAt(i,j); T sum = this->ElementAt(i,j);
for(k = 0; k < j; k++) for(k = 0; k < j; k++)
sum -= this->ElementAt(i,k)*this->ElementAt(k,j); sum -= this->ElementAt(i,k)*this->ElementAt(k,j);
this->ElementAt(i,j) = sum; this->ElementAt(i,j) = sum;
T t = scaling[i] * math::Abs(sum); T t = scaling[i] * math::Abs(sum);
if(t >= largest) { if(t >= largest) {
largest = t; largest = t;
imax = i; imax = i;
} }
} }
if (j != imax) { if (j != imax) {
for (int k = 0; k < 4; k++) { for (int k = 0; k < 4; k++) {
T dum = this->ElementAt(imax,k); T dum = this->ElementAt(imax,k);
this->ElementAt(imax,k) = this->ElementAt(j,k); this->ElementAt(imax,k) = this->ElementAt(j,k);
this->ElementAt(j,k) = dum; this->ElementAt(j,k) = dum;
} }
d = -(d); d = -(d);
scaling[imax] = scaling[j]; scaling[imax] = scaling[j];
} }
index[j]=imax; index[j]=imax;
if (this->ElementAt(j,j) == 0.0) this->ElementAt(j,j) = (T)TINY; if (this->ElementAt(j,j) == 0.0) this->ElementAt(j,j) = (T)TINY;
if (j != 3) { if (j != 3) {
T dum = (T)1.0 / (this->ElementAt(j,j)); T dum = (T)1.0 / (this->ElementAt(j,j));
for (i = j+1; i < 4; i++) for (i = j+1; i < 4; i++)
this->ElementAt(i,j) *= dum; this->ElementAt(i,j) *= dum;
} }
} }
return true; return true;
} }
template <class T> Point4<T> LinearSolve<T>::Solve(const Point4<T> &b) { template <class T> Point4<T> LinearSolve<T>::Solve(const Point4<T> &b) {
Point4<T> x(b); Point4<T> x(b);
int first = -1, ip; int first = -1, ip;
for(int i = 0; i < 4; i++) { for(int i = 0; i < 4; i++) {
ip = index[i]; ip = index[i];
T sum = x[ip]; T sum = x[ip];
x[ip] = x[i]; x[ip] = x[i];
if(first!= -1) if(first!= -1)
for(int j = first; j <= i-1; j++) for(int j = first; j <= i-1; j++)
sum -= this->ElementAt(i,j) * x[j]; sum -= this->ElementAt(i,j) * x[j];
else else
if(sum) first = i; if(sum) first = i;
x[i] = sum; x[i] = sum;
} }
for (int i = 3; i >= 0; i--) { for (int i = 3; i >= 0; i--) {
T sum = x[i]; T sum = x[i];
for (int j = i+1; j < 4; j++) for (int j = i+1; j < 4; j++)
sum -= this->ElementAt(i, j) * x[j]; sum -= this->ElementAt(i, j) * x[j];
x[i] = sum / this->ElementAt(i, i); x[i] = sum / this->ElementAt(i, i);
} }
return x; return x;
} }

57
vcg/math/eigen.h
View File

@ -31,7 +31,24 @@
// forward declarations // forward declarations
namespace Eigen { namespace Eigen {
#include "../Eigen/src/Core/util/Meta.h"
template<typename Derived1, typename Derived2, int Size> struct ei_lexi_comparison; template<typename Derived1, typename Derived2, int Size> struct ei_lexi_comparison;
template<typename Derived1, typename Derived2,
bool SameType = ei_is_same_type<Derived1,Derived2>::ret,
bool SameSize = Derived1::SizeAtCompileTime==Derived2::SizeAtCompileTime>
struct ei_import_selector;
template<typename XprType,
int Rows = XprType::RowsAtCompileTime,
int Cols = XprType::ColsAtCompileTime,
int StorageOrder = XprType::Flags&1,
int MRows = XprType::MaxRowsAtCompileTime,
int MCols = XprType::MaxColsAtCompileTime>
struct ei_to_vcgtype;
} }
#include "base.h" #include "base.h"
@ -156,6 +173,46 @@ template<typename Derived1, typename Derived2> struct ei_lexi_comparison<Derived
} }
}; };
// implementation of Import
template<typename Derived1, typename Derived2>
struct ei_import_selector<Derived1,Derived2,true,true>
{
static void run(Derived1& a, const Derived2& b) { a = b; }
};
template<typename Derived1, typename Derived2>
struct ei_import_selector<Derived1,Derived2,false,true>
{
static void run(Derived1& a, const Derived2& b)
{ a = b.template cast<typename Derived1::Scalar>(); }
};
template<typename Derived1, typename Derived2>
struct ei_import_selector<Derived1,Derived2,false,false>
{
static void run(Derived1& a, const Derived2& b)
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived1);
EIGEN_STATIC_ASSERT_FIXED_SIZE(Derived1);
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived2);
EIGEN_STATIC_ASSERT_FIXED_SIZE(Derived2);
enum {
Size1 = Derived1::SizeAtCompileTime,
Size2 = Derived2::SizeAtCompileTime
};
assert(Size1<=4 && Size2<=4);
a.coeffRef(0) = Scalar(b.coeff(0));
if (Size1>1) { if (Size2>1) a.coeffRef(1) = Scalar(b.coeff(1)); else a.coeffRef(1) = 0; }
if (Size1>2) { if (Size2>2) a.coeffRef(2) = Scalar(b.coeff(2)); else a.coeffRef(2) = 0; }
if (Size1>3) { if (Size2>3) a.coeffRef(3) = Scalar(b.coeff(3)); else a.coeffRef(3) = 0; }
}
};
// default implementation of ei_to_vcgtype
// the specialization are with
template<typename XprType,int Rows,int Cols,int StorageOrder,int MRows,int MCols>
struct ei_to_vcgtype { typedef Matrix<typename XprType::Scalar,Rows,Cols,StorageOrder,MRows,MCols> type; };
} }
#define VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATOR(Derived, Op) \ #define VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATOR(Derived, Op) \

54
vcg/math/eigen_matrix_addons.h
View File

@ -24,28 +24,31 @@
#warning You are including deprecated math stuff #warning You are including deprecated math stuff
enum {Dimension = SizeAtCompileTime}; enum {Dimension = SizeAtCompileTime};
typedef typename ei_to_vcgtype<Matrix>::type EquivVcgType;
typedef vcg::VoidType ParamType; typedef vcg::VoidType ParamType;
typedef Matrix PointType; typedef Matrix PointType;
using Base::V; using Base::V;
// automatic conversion to similar vcg types
// the otherway round is implicit because they inherits this Matrix tyoe
operator EquivVcgType& () { return *reinterpret_cast<EquivVcgType*>(this); }
operator const EquivVcgType& () const { return *reinterpret_cast<const EquivVcgType*>(this); }
/** \deprecated use m.cast<NewScalar>() */
/// importer for points with different scalar type and-or dimensionality /// importer for points with different scalar type and-or dimensionality
// FIXME the Point3/Point4 specialization were only for same sizes ?? // FIXME the Point3/Point4 specialization were only for same sizes ??
// while the Point version was generic like this one // while the Point version was generic like this one
template<typename OtherDerived> template<typename OtherDerived>
inline void Import(const MatrixBase<OtherDerived>& b) inline void Import(const MatrixBase<OtherDerived>& b)
{ {
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix); ei_import_selector<Matrix,OtherDerived>::run(*this,b.derived());
EIGEN_STATIC_ASSERT_FIXED_SIZE(Matrix);
EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived);
EIGEN_STATIC_ASSERT_FIXED_SIZE(OtherDerived);
enum { OtherSize = OtherDerived::SizeAtCompileTime };
assert(SizeAtCompileTime<=4 && OtherSize<=4);
data()[0] = Scalar(b[0]);
if (SizeAtCompileTime>1) { if (OtherSize>1) data()[1] = Scalar(b[1]); else data()[1] = 0; }
if (SizeAtCompileTime>2) { if (OtherSize>2) data()[2] = Scalar(b[2]); else data()[2] = 0; }
if (SizeAtCompileTime>3) { if (OtherSize>3) data()[3] = Scalar(b[3]); else data()[3] = 0; }
} }
/// constructor for points with different scalar type and-or dimensionality
template<typename OtherDerived>
static inline Matrix Construct(const MatrixBase<OtherDerived>& b)
{ Matrix p; p.Import(b); return p; }
/// importer for homogeneous points /// importer for homogeneous points
template<typename OtherDerived> template<typename OtherDerived>
inline void ImportHomo(const MatrixBase<OtherDerived>& b) inline void ImportHomo(const MatrixBase<OtherDerived>& b)
@ -53,26 +56,15 @@ inline void ImportHomo(const MatrixBase<OtherDerived>& b)
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix); EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix);
EIGEN_STATIC_ASSERT_FIXED_SIZE(Matrix); EIGEN_STATIC_ASSERT_FIXED_SIZE(Matrix);
EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,SizeAtCompileTime-1); EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,SizeAtCompileTime-1);
this->template start<SizeAtCompileTime-1> = b; this->template start<SizeAtCompileTime-1> = b;
data()[SizeAtCompileTime-1] = Scalar(1.0); data()[SizeAtCompileTime-1] = Scalar(1.0);
} }
/// constructor for points with different scalar type and-or dimensionality /// constructor for homogeneus point.
template<typename OtherDerived>
static inline Matrix Construct(const MatrixBase<OtherDerived>& b)
{
Matrix p; p.Import(b);
return p;
}
/// constructor for homogeneus point.
template<typename OtherDerived> template<typename OtherDerived>
static inline Matrix ConstructHomo(const MatrixBase<OtherDerived>& b) static inline Matrix ConstructHomo(const MatrixBase<OtherDerived>& b)
{ { Matrix p; p.ImportHomo(b); return p; }
Matrix p; p.ImportHomo(b);
return p;
}
inline const Scalar &X() const { return data()[0]; } inline const Scalar &X() const { return data()[0]; }
inline const Scalar &Y() const { return data()[1]; } inline const Scalar &Y() const { return data()[1]; }
@ -81,15 +73,19 @@ inline Scalar &X() { return data()[0]; }
inline Scalar &Y() { return data()[1]; } inline Scalar &Y() { return data()[1]; }
inline Scalar &Z() { assert(SizeAtCompileTime>2); return data()[2]; } inline Scalar &Z() { assert(SizeAtCompileTime>2); return data()[2]; }
// note, W always returns the last entry /** note, W always returns the last entry */
inline Scalar& W() { return data()[SizeAtCompileTime-1]; } inline Scalar& W() { return data()[SizeAtCompileTime-1]; }
/** note, W always returns the last entry */
inline const Scalar& W() const { return data()[SizeAtCompileTime-1]; } inline const Scalar& W() const { return data()[SizeAtCompileTime-1]; }
Scalar* V() { return data(); } /** \deprecated use .data() */
const Scalar* V() const { return data(); } EIGEN_DEPRECATED Scalar* V() { return data(); }
/** \deprecated use .data() */
EIGEN_DEPRECATED const Scalar* V() const { return data(); }
/** \deprecated use m.coeff(i) or m[i] or m(i) */
// overloaded to return a const reference // overloaded to return a const reference
inline const Scalar& V( const int i ) const EIGEN_DEPRECATED inline const Scalar& V( const int i ) const
{ {
assert(i>=0 && i<SizeAtCompileTime); assert(i>=0 && i<SizeAtCompileTime);
return data()[i]; return data()[i];
@ -106,7 +102,7 @@ inline Matrix LocalToGlobal(ParamType p) const { return *this; }
* (provided for uniformity with other spatial classes. trivial for points) */ * (provided for uniformity with other spatial classes. trivial for points) */
inline ParamType GlobalToLocal(PointType /*p*/) const { return ParamType(); } inline ParamType GlobalToLocal(PointType /*p*/) const { return ParamType(); }
/** /**
* Convert to polar coordinates from cartesian coordinates. * Convert to polar coordinates from cartesian coordinates.
* *
* Theta is the azimuth angle and ranges between [0, 360) degrees. * Theta is the azimuth angle and ranges between [0, 360) degrees.

5
vcg/math/matrix.h
View File

@ -41,6 +41,8 @@ template<class Scalar> class Matrix;
namespace Eigen{ namespace Eigen{
template<typename Scalar> template<typename Scalar>
struct ei_traits<vcg::ndim::Matrix<Scalar> > : ei_traits<Eigen::Matrix<Scalar,Dynamic,Dynamic> > {}; struct ei_traits<vcg::ndim::Matrix<Scalar> > : ei_traits<Eigen::Matrix<Scalar,Dynamic,Dynamic> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,Dynamic,Dynamic,RowMajor,Dynamic,Dynamic>
{ typedef vcg::ndim::Matrix<typename XprType::Scalar> type; };
} }
namespace vcg{ namespace vcg{
@ -50,7 +52,7 @@ namespace ndim{
/* @{ */ /* @{ */
/*! /*!
* \deprecated use Matrix<Scalar,Rows,Cols> or Matrix<Scalar,Dynamic,Dynamic> * \deprecated use Matrix<Scalar,Rows,Cols> or Matrix<Scalar,Dynamic,Dynamic> or any typedef
* This class represent a generic <I>m</I><EFBFBD><I>n</I> matrix. The class is templated over the scalar type field. * This class represent a generic <I>m</I><EFBFBD><I>n</I> matrix. The class is templated over the scalar type field.
* @param Scalar (Templete Parameter) Specifies the ScalarType field. * @param Scalar (Templete Parameter) Specifies the ScalarType field.
*/ */
@ -63,7 +65,6 @@ public:
_EIGEN_GENERIC_PUBLIC_INTERFACE(Matrix,_Base); _EIGEN_GENERIC_PUBLIC_INTERFACE(Matrix,_Base);
typedef _Scalar ScalarType; typedef _Scalar ScalarType;
VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Matrix) VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Matrix)
/*! /*!

20
vcg/math/matrix33.h
View File

@ -38,12 +38,14 @@ template<class Scalar> class Matrix33;
namespace Eigen{ namespace Eigen{
template<typename Scalar> template<typename Scalar>
struct ei_traits<vcg::Matrix33<Scalar> > : ei_traits<Eigen::Matrix<Scalar,3,3,RowMajor> > {}; struct ei_traits<vcg::Matrix33<Scalar> > : ei_traits<Eigen::Matrix<Scalar,3,3,RowMajor> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,3,3,RowMajor,3,3>
{ typedef vcg::Matrix33<typename XprType::Scalar> type; };
} }
namespace vcg { namespace vcg {
/** \deprecated use Matrix<Scalar,3,3>
/** @name Matrix33 @name Matrix33
Class Matrix33. Class Matrix33.
This is the class for definition of a matrix 3x3. This is the class for definition of a matrix 3x3.
@param S (Templete Parameter) Specifies the ScalarType field. @param S (Templete Parameter) Specifies the ScalarType field.
@ -98,7 +100,6 @@ public:
return SetRotateRad(math::ToRad(angle),axis); return SetRotateRad(math::ToRad(angle),axis);
} }
// Warning, this Inversion code can be HIGHLY NUMERICALLY UNSTABLE! // Warning, this Inversion code can be HIGHLY NUMERICALLY UNSTABLE!
// In most case you are advised to use the Invert() method based on SVD decomposition. // In most case you are advised to use the Invert() method based on SVD decomposition.
/** \deprecated */ /** \deprecated */
@ -110,18 +111,13 @@ public:
printf("| %g \t%g \t%g |\n",coeff(i,0),coeff(i,1),coeff(i,2)); printf("| %g \t%g \t%g |\n",coeff(i,0),coeff(i,1),coeff(i,2));
} }
/* /** \deprecated use a * b.transpose()
compute the matrix generated by the product of a * b^T compute the matrix generated by the product of a * b^T
*/ */
void ExternalProduct(const Point3<Scalar> &a, const Point3<Scalar> &b) // hm.... this is the outer product
{ void ExternalProduct(const Point3<Scalar> &a, const Point3<Scalar> &b) { *this = a * b.transpose(); }
for(int i=0;i<3;++i)
for(int j=0;j<3;++j)
(*this)[i][j] = a[i]*b[j];
}
/* Compute the Frobenius Norm of the Matrix /** Compute the Frobenius Norm of the Matrix */
*/
Scalar Norm() { return Base::cwise().abs2().sum(); } Scalar Norm() { return Base::cwise().abs2().sum(); }
// { // {
// // FIXME looks like there was a bug: j is not used !!! // // FIXME looks like there was a bug: j is not used !!!

65
vcg/math/matrix44.h
View File

@ -41,6 +41,8 @@ template<class Scalar> class Matrix44;
namespace Eigen{ namespace Eigen{
template<typename Scalar> template<typename Scalar>
struct ei_traits<vcg::Matrix44<Scalar> > : ei_traits<Eigen::Matrix<Scalar,4,4,RowMajor> > {}; struct ei_traits<vcg::Matrix44<Scalar> > : ei_traits<Eigen::Matrix<Scalar,4,4,RowMajor> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,4,4,RowMajor,4,4>
{ typedef vcg::Matrix44<typename XprType::Scalar> type; };
} }
namespace vcg { namespace vcg {
@ -140,17 +142,6 @@ public:
Matrix44 &SetRotateDeg(Scalar AngleDeg, const Point3<Scalar> & axis); Matrix44 &SetRotateDeg(Scalar AngleDeg, const Point3<Scalar> & axis);
Matrix44 &SetRotateRad(Scalar AngleRad, const Point3<Scalar> & axis); Matrix44 &SetRotateRad(Scalar AngleRad, const Point3<Scalar> & axis);
template <class Q> void Import(const Matrix44<Q> &m) {
for(int i = 0; i < 16; i++)
Base::data()[i] = (Scalar)(m.data()[i]);
}
template <class Q>
static inline Matrix44 Construct( const Matrix44<Q> & b )
{
Matrix44 tmp; tmp.FromMatrix(b);
return tmp;
}
// template <class T> Point3<T> operator*(const Point3<T> &p) { // template <class T> Point3<T> operator*(const Point3<T> &p) {
// T w; // T w;
// Point3<T> s; // Point3<T> s;
@ -173,30 +164,10 @@ public:
return s; return s;
} }
void print() {std::cout << *this << "\n\n";}
}; };
/** Class for solving A * x = b. */
template <class T> class LinearSolve: public Matrix44<T> {
public:
LinearSolve(const Matrix44<T> &m);
Point4<T> Solve(const Point4<T> &b); // solve A <20> x = b
///If you need to solve some equation you can use this function instead of Matrix44 one for speed.
T Determinant() const;
protected:
///Holds row permutation.
int index[4]; //hold permutation
///Hold sign of row permutation (used for determinant sign)
T d;
bool Decompose();
};
/*** Postmultiply */
//template <class T> Point3<T> operator*(const Point3<T> &p, const Matrix44<T> &m);
///Premultiply
// template <class T> Point3<T> operator*(const Matrix44<T> &m, const Point3<T> &p);
//return NULL matrix if not invertible //return NULL matrix if not invertible
template <class T> Matrix44<T> &Invert(Matrix44<T> &m); template <class T> Matrix44<T> &Invert(Matrix44<T> &m);
template <class T> Matrix44<T> Inverse(const Matrix44<T> &m); template <class T> Matrix44<T> Inverse(const Matrix44<T> &m);
@ -206,18 +177,6 @@ typedef Matrix44<int> Matrix44i;
typedef Matrix44<float> Matrix44f; typedef Matrix44<float> Matrix44f;
typedef Matrix44<double> Matrix44d; typedef Matrix44<double> Matrix44d;
//template <class T> T &Matrix44<T>::operator[](const int i) {
// assert(i >= 0 && i < 16);
// return ((T *)_a)[i];
//}
//
//template <class T> const T &Matrix44<T>::operator[](const int i) const {
// assert(i >= 0 && i < 16);
// return ((T *)_a)[i];
//}
template < class PointType , class T > void operator*=( std::vector<PointType> &vert, const Matrix44<T> & m ) { template < class PointType , class T > void operator*=( std::vector<PointType> &vert, const Matrix44<T> & m ) {
typename std::vector<PointType>::iterator ii; typename std::vector<PointType>::iterator ii;
for(ii=vert.begin();ii!=vert.end();++ii) for(ii=vert.begin();ii!=vert.end();++ii)
@ -487,22 +446,6 @@ bool Decompose(Matrix44<T> &M, Point3<T> &ScaleV, Point3<T> &ShearV, Point3<T> &
return true; return true;
} }
// template <class T> Point3<T> operator*(const Point3<T> &p, const Matrix44<T> &m) {
// T w;
// Point3<T> s;
// s[0] = m.ElementAt(0, 0)*p[0] + m.ElementAt(1, 0)*p[1] + m.ElementAt(2, 0)*p[2] + m.ElementAt(3, 0);
// s[1] = m.ElementAt(0, 1)*p[0] + m.ElementAt(1, 1)*p[1] + m.ElementAt(2, 1)*p[2] + m.ElementAt(3, 1);
// s[2] = m.ElementAt(0, 2)*p[0] + m.ElementAt(1, 2)*p[1] + m.ElementAt(2, 2)*p[2] + m.ElementAt(3, 2);
// w = m.ElementAt(0, 3)*p[0] + m.ElementAt(1, 3)*p[1] + m.ElementAt(2, 3)*p[2] + m.ElementAt(3, 3);
// if(w != 0) s /= w;
// return s;
// }
/* /*
To invert a matrix you can To invert a matrix you can
either invert the matrix inplace calling either invert the matrix inplace calling

40
vcg/math/shot.h
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -93,9 +93,9 @@ creation
****************************************************************************/ ****************************************************************************/
/** class Shot /** class Shot
The shot is made of two things: The shot is made of two things:
* the Instrinsics paramaters, which are stored as a Camera type (check vcg/math/camera) and that * the Instrinsics paramaters, which are stored as a Camera type (check vcg/math/camera) and that
determines how a point in the frame of the camera is projected in the 2D projection plane determines how a point in the frame of the camera is projected in the 2D projection plane
* the Extrinsics parameters, which are stored in the class Shot (che the type ReferenceFrame) * the Extrinsics parameters, which are stored in the class Shot (che the type ReferenceFrame)
@ -103,7 +103,7 @@ and that tell viewpoint and view direction.
The Extrinsics parameters are kept as a rotation matrix "rot" and a translation vector "tra" The Extrinsics parameters are kept as a rotation matrix "rot" and a translation vector "tra"
NOTE: the translation matrix "tra" corresponds to -viewpoint while the rotation matrix NOTE: the translation matrix "tra" corresponds to -viewpoint while the rotation matrix
"rot" corresponds to the axis of the reference frame by row, i.e. "rot" corresponds to the axis of the reference frame by row, i.e.
rot[0][0|1|2] == X axis rot[0][0|1|2] == X axis
rot[1][0|1|2] == Y axis rot[1][0|1|2] == Y axis
rot[2][0|1|2] == Z axis rot[2][0|1|2] == Z axis
@ -136,7 +136,7 @@ public:
class ReferenceFrame { class ReferenceFrame {
friend class Shot<ScalarType, RotoType>; friend class Shot<ScalarType, RotoType>;
RotoType rot; // rotation RotoType rot; // rotation
Point3<S> tra; // viewpoint Point3<S> tra; // viewpoint
public: public:
void SetIdentity(){ rot.SetIdentity(); tra = Point3<S>(0.0,0.0,0.0);} void SetIdentity(){ rot.SetIdentity(); tra = Point3<S>(0.0,0.0,0.0);}
void SetTra(const Point3<S> & tr) {tra = tr;} void SetTra(const Point3<S> & tr) {tra = tr;}
@ -146,7 +146,7 @@ public:
}; };
Camera<S> Intrinsics; // the camera that made the shot Camera<S> Intrinsics; // the camera that made the shot
ReferenceFrame<S,RotationType> Extrinsics; // the position and orientation of the camera ReferenceFrame<S,RotationType> Extrinsics; // the position and orientation of the camera
Shot(Camera<S> c) Shot(Camera<S> c)
@ -186,13 +186,13 @@ public:
/// convert a 3d point from camera to world coordinates /// convert a 3d point from camera to world coordinates
vcg::Point3<S> ConvertCameraToWorldCoordinates(const vcg::Point3<S> & p) const; vcg::Point3<S> ConvertCameraToWorldCoordinates(const vcg::Point3<S> & p) const;
/// project a 3d point from world coordinates to 2d camera viewport (pixel) /// project a 3d point from world coordinates to 2d camera viewport (pixel)
vcg::Point2<S> Project(const vcg::Point3<S> & p) const; vcg::Point2<S> Project(const vcg::Point3<S> & p) const;
/// inverse projection from 2d camera viewport (pixel) + Zdepth to 3d world coordinates /// inverse projection from 2d camera viewport (pixel) + Zdepth to 3d world coordinates
vcg::Point3<S> UnProject(const vcg::Point2<S> & p, const S & d) const; vcg::Point3<S> UnProject(const vcg::Point2<S> & p, const S & d) const;
/// returns distance of point p from camera plane (z depth), used for unprojection /// returns distance of point p from camera plane (z depth), used for unprojection
S Depth(const vcg::Point3<S> & p)const; S Depth(const vcg::Point3<S> & p)const;
@ -205,8 +205,7 @@ public:
Matrix44<S> GetExtrinsicsToWorldMatrix() const { Matrix44<S> GetExtrinsicsToWorldMatrix() const {
Matrix44<S> rotM ; Matrix44<S> rotM ;
Extrinsics.rot.ToMatrix(rotM); Extrinsics.rot.ToMatrix(rotM);
Transpose(rotM); return Matrix44<S>().SetTranslate(Extrinsics.tra) * rotM.transpose();
return Matrix44<S>().SetTranslate(Extrinsics.tra) * rotM;
} }
/* Returns the matrix M such that /* Returns the matrix M such that
@ -240,7 +239,7 @@ public:
//--- //---
/// GET the viewpoint /// GET the viewpoint
template <class S, class RotationType> template <class S, class RotationType>
const vcg::Point3<S> Shot<S,RotationType>::GetViewPoint() const const vcg::Point3<S> Shot<S,RotationType>::GetViewPoint() const
{ {
return Extrinsics.tra; return Extrinsics.tra;
} }
@ -254,10 +253,10 @@ void Shot<S,RotationType>::SetViewPoint(const vcg::Point3<S> & viewpoint)
/// GET the i-th axis of the coordinate system of the camera /// GET the i-th axis of the coordinate system of the camera
template <class S, class RotationType> template <class S, class RotationType>
vcg::Point3<S> Shot<S,RotationType>::Axis(const int & i) const vcg::Point3<S> Shot<S,RotationType>::Axis(const int & i) const
{ {
vcg::Matrix44<S> m; vcg::Matrix44<S> m;
Extrinsics.rot.ToMatrix(m); Extrinsics.rot.ToMatrix(m);
vcg::Point3<S> aa = m.GetRow3(i); vcg::Point3<S> aa = m.GetRow3(i);
return aa; return aa;
} }
@ -285,7 +284,7 @@ void Shot<S,RotationType>::LookTowards(const vcg::Point3<S> & z_dir,const vcg::P
{ {
vcg::Point3<S> x_dir = up ^-z_dir; vcg::Point3<S> x_dir = up ^-z_dir;
vcg::Point3<S> y_dir = -z_dir ^x_dir; vcg::Point3<S> y_dir = -z_dir ^x_dir;
Matrix44<S> m; Matrix44<S> m;
m.SetIdentity(); m.SetIdentity();
*(vcg::Point3<S> *)&m[0][0] = x_dir/x_dir.Norm(); *(vcg::Point3<S> *)&m[0][0] = x_dir/x_dir.Norm();
@ -316,12 +315,11 @@ vcg::Point3<S> Shot<S,RotationType>::ConvertCameraToWorldCoordinates(const vcg::
vcg::Point3<S> cp = p; vcg::Point3<S> cp = p;
cp[2]=-cp[2]; // note: World reference system is left handed cp[2]=-cp[2]; // note: World reference system is left handed
Extrinsics.rot.ToMatrix(rotM); Extrinsics.rot.ToMatrix(rotM);
Transpose(rotM); cp = rotM.transpose() * cp + GetViewPoint();
cp = rotM * cp+ GetViewPoint();
return cp; return cp;
} }
/// project a 3d point from world coordinates to 2d camera viewport (pixel) /// project a 3d point from world coordinates to 2d camera viewport (pixel)
template <class S, class RotationType> template <class S, class RotationType>
vcg::Point2<S> Shot<S,RotationType>::Project(const vcg::Point3<S> & p) const vcg::Point2<S> Shot<S,RotationType>::Project(const vcg::Point3<S> & p) const
{ {
@ -341,9 +339,9 @@ vcg::Point3<S> Shot<S,RotationType>::UnProject(const vcg::Point2<S> & p, const S
return wp; return wp;
} }
/// returns distance of point p from camera plane (z depth), used for unprojection /// returns distance of point p from camera plane (z depth), used for unprojection
template <class S, class RotationType> template <class S, class RotationType>
S Shot<S,RotationType>::Depth(const vcg::Point3<S> & p)const S Shot<S,RotationType>::Depth(const vcg::Point3<S> & p)const
{ {
return ConvertWorldToCameraCoordinates(p).Z(); return ConvertWorldToCameraCoordinates(p).Z();
} }

2
vcg/space/point2.h
View File

@ -37,6 +37,8 @@ template<typename Scalar> class Point2;
namespace Eigen { namespace Eigen {
template<typename Scalar> struct ei_traits<vcg::Point2<Scalar> > : ei_traits<Eigen::Matrix<Scalar,2,1> > {}; template<typename Scalar> struct ei_traits<vcg::Point2<Scalar> > : ei_traits<Eigen::Matrix<Scalar,2,1> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,2,1,0,2,1>
{ typedef vcg::Point2<typename XprType::Scalar> type; };
} }
namespace vcg { namespace vcg {

9
vcg/space/point3.h
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -38,6 +38,9 @@ namespace Eigen{
template<typename Scalar> struct ei_traits<vcg::Point3<Scalar> > : ei_traits<Eigen::Matrix<Scalar,3,1> > {}; template<typename Scalar> struct ei_traits<vcg::Point3<Scalar> > : ei_traits<Eigen::Matrix<Scalar,3,1> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,3,1,0,3,1>
{ typedef vcg::Point3<typename XprType::Scalar> type; };
template<typename Scalar> template<typename Scalar>
struct NumTraits<vcg::Point3<Scalar> > : NumTraits<Scalar> struct NumTraits<vcg::Point3<Scalar> > : NumTraits<Scalar>
{ {
@ -79,7 +82,7 @@ public:
_EIGEN_GENERIC_PUBLIC_INTERFACE(Point3,_Base); _EIGEN_GENERIC_PUBLIC_INTERFACE(Point3,_Base);
VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Point3) VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Point3)
/** @name Standard Constructors and Initializers /** @name Standard Constructors and Initializers
No casting operators have been introduced to avoid automatic unattended (and costly) conversion between different point types No casting operators have been introduced to avoid automatic unattended (and costly) conversion between different point types
**/ **/
@ -92,7 +95,7 @@ public:
// this one is very useless // this one is very useless
template <class Q> template <class Q>
static inline Point3 Construct( const Q & P0, const Q & P1, const Q & P2) static inline Point3 Construct( const Q & P0, const Q & P1, const Q & P2)
{ {
return Point3(Scalar(P0),Scalar(P1),Scalar(P2)); return Point3(Scalar(P0),Scalar(P1),Scalar(P2));

4
vcg/space/point4.h
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -36,6 +36,8 @@ template<typename Scalar> class Point4;
namespace Eigen { namespace Eigen {
template<typename Scalar> struct ei_traits<vcg::Point4<Scalar> > : ei_traits<Eigen::Matrix<Scalar,4,1> > {}; template<typename Scalar> struct ei_traits<vcg::Point4<Scalar> > : ei_traits<Eigen::Matrix<Scalar,4,1> > {};
template<typename XprType> struct ei_to_vcgtype<XprType,4,1,0,4,1>
{ typedef vcg::Point4<typename XprType::Scalar> type; };
} }
namespace vcg { namespace vcg {

144
wrap/gl/deprecated_math.h Normal file
View File

@ -0,0 +1,144 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2004 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
/****************************************************************************
History
$Log: not supported by cvs2svn $
Revision 1.10 2006/02/13 13:05:05 cignoni
Removed glew inclusion
Revision 1.9 2004/09/30 00:48:07 ponchio
<gl/glew.h> -> <GL/glew.h>
Revision 1.8 2004/09/28 14:04:36 ganovelli
glGet added
Revision 1.7 2004/07/13 15:55:57 cignoni
Added test on presence of glTranspose extension (for old hw support)
Revision 1.6 2004/05/26 15:12:39 cignoni
Removed inclusion of gl extension stuff
Revision 1.5 2004/05/12 20:54:55 ponchio
*** empty log message ***
Revision 1.4 2004/05/12 13:07:47 ponchio
Added #include <glew.h>
Revision 1.3 2004/05/04 23:36:23 cignoni
remove include of gl and added glextgension exploiting,
Revision 1.2 2004/04/07 10:47:03 cignoni
inlined functions for avoid multiple linking errors
Revision 1.1 2004/03/31 15:27:17 ponchio
*** empty log message ***
****************************************************************************/
#ifndef VCG_GL_MATH_H
#define VCG_GL_MATH_H
// Please note that this file assume that you have already included your
// gl-extension wrapping utility, and that therefore all the extension symbol are already defined.
#include <vcg/math/matrix44.h>
#include <vcg/math/similarity.h>
//#include <GL/glew.h> // please do not include it!
namespace vcg {
inline void glMultMatrixE(const Matrix44f &matrix) {
//glMultMatrixf((const GLfloat *)(matrix[0]));
if(glMultTransposeMatrixf) glMultTransposeMatrixf((const GLfloat *)(matrix.V()));
else {
glMultMatrixf((const GLfloat *)(matrix.transpose().eval().V()));
}
}
inline void glMultMatrix(const Matrix44f &matrix) {
glMultMatrixf((const GLfloat *)(matrix.transpose().eval().V()));
}
inline void glMultMatrixE(const Matrix44d &matrix) {
if(glMultTransposeMatrixd) glMultTransposeMatrixd((const GLdouble *)(matrix.V()));
else {
glMultMatrixd((const GLdouble *)(matrix.transpose().eval().V()));
}
}
inline void glMultMatrix(const Matrix44d &matrix) {
glMultMatrixd((const GLdouble *)(matrix.transpose().eval().V()));
}
inline void glMultMatrixDirect(const Matrix44f &matrix) {
glMultMatrixf((const GLfloat *)(matrix.V()));
}
inline void glMultMatrixDirect(const Matrix44d &matrix) {
glMultMatrixd((const GLdouble *)(matrix.V()));
}
inline void glMultMatrix(const Similarityf &s) {
glTranslatef(s.tra[0], s.tra[1], s.tra[2]);
glScalef(s.sca, s.sca, s.sca);
float alpha;
Point3f axis;
s.rot.ToAxis(alpha, axis);
glRotatef(math::ToDeg(alpha), axis[0], axis[1], axis[2]);
}
inline void glMultMatrix(const Similarityd &s) {
glTranslated(s.tra[0], s.tra[1], s.tra[2]);
double alpha;
Point3d axis;
s.rot.ToAxis(alpha, axis);
glRotated(math::ToDeg(alpha), axis[0], axis[1], axis[2]);
glScaled(s.sca, s.sca, s.sca);
}
inline void glGetv(const GLenum pname, Matrix44f & m){
Matrix44f tmp;
glGetFloatv(pname,tmp.V());
m = tmp.transpose();
}
inline void glGetv(const GLenum pname, Matrix44d & m){
Matrix44d tmp;
glGetDoublev(pname,tmp.V());
m = tmp.transpose();
}
inline void glGetDirectv(const GLenum pname, Matrix44f & m){
glGetFloatv(pname,m.V());
}
inline void glGetDirecv(const GLenum pname, Matrix44d & m){
glGetDoublev(pname,m.V());
}
}//namespace
#endif

90
wrap/gl/math.h
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -57,10 +57,14 @@ Revision 1.1 2004/03/31 15:27:17 ponchio
****************************************************************************/ ****************************************************************************/
#ifndef VCG_USE_EIGEN
#include "deprecated_math.h"
#else
#ifndef VCG_GL_MATH_H #ifndef VCG_GL_MATH_H
#define VCG_GL_MATH_H #define VCG_GL_MATH_H
// Please note that this file assume that you have already included your // Please note that this file assume that you have already included your
// gl-extension wrapping utility, and that therefore all the extension symbol are already defined. // gl-extension wrapping utility, and that therefore all the extension symbol are already defined.
#include <vcg/math/matrix44.h> #include <vcg/math/matrix44.h>
@ -70,44 +74,41 @@ Revision 1.1 2004/03/31 15:27:17 ponchio
namespace vcg { namespace vcg {
inline void glMultMatrixE(const Matrix44f &matrix) { inline void glLoadMatrix(const Eigen::Matrix<float,4,4>& matrix) { glLoadMatrixf(matrix.data()); }
//glMultMatrixf((const GLfloat *)(matrix[0])); inline void glLoadMatrix(const Eigen::Matrix<float,4,4,Eigen::RowMajor>& matrix) {
if(glMultTransposeMatrixf) glMultTransposeMatrixf((const GLfloat *)(matrix.V())); Eigen::Matrix4f tmp(matrix);
else { glLoadMatrixf(tmp.data());
glMultMatrixf((const GLfloat *)(matrix.transpose().eval().V()));
}
} }
inline void glLoadMatrix(const Eigen::Matrix<double,4,4>& matrix) { glLoadMatrixd(matrix.data()); }
inline void glLoadMatrix(const Eigen::Matrix<double,4,4,Eigen::RowMajor>& matrix) {
Eigen::Matrix4d tmp(matrix);
glLoadMatrixd(tmp.data());
}
template<typename Scalar>
inline void glLoadMatrix(const Eigen::Transform<Scalar,3>& t) { glLoadMatrix(t.matrix()); }
inline void glMultMatrix(const Matrix44f &matrix) {
glMultMatrixf((const GLfloat *)(matrix.transpose().eval().V()));
}
inline void glMultMatrixE(const Matrix44d &matrix) { inline void glMultMatrix(const Eigen::Matrix<float,4,4>& matrix) { glMultMatrixf(matrix.data()); }
if(glMultTransposeMatrixd) glMultTransposeMatrixd((const GLdouble *)(matrix.V())); inline void glMultMatrix(const Eigen::Matrix<float,4,4,Eigen::RowMajor>& matrix) {
else { Eigen::Matrix<float,4,4> tmp(matrix);
glMultMatrixd((const GLdouble *)(matrix.transpose().eval().V())); glMultMatrixf(tmp.data());
}
} }
inline void glMultMatrix(const Matrix44d &matrix) { inline void glMultMatrix(const Eigen::Matrix<double,4,4>& matrix) { glMultMatrixd(matrix.data()); }
glMultMatrixd((const GLdouble *)(matrix.transpose().eval().V())); inline void glMultMatrix(const Eigen::Matrix<double,4,4,Eigen::RowMajor>& matrix) {
Eigen::Matrix<double,4,4> tmp(matrix);
glMultMatrixd(tmp.data());
} }
template<typename Scalar>
inline void glMultMatrix(const Eigen::Transform<Scalar,3>& t) { glMultMatrix(t.matrix()); }
inline void glMultMatrixDirect(const Matrix44f &matrix) {
glMultMatrixf((const GLfloat *)(matrix.V()));
}
inline void glMultMatrixDirect(const Matrix44d &matrix) {
glMultMatrixd((const GLdouble *)(matrix.V()));
}
inline void glMultMatrix(const Similarityf &s) { inline void glMultMatrix(const Similarityf &s) {
glTranslatef(s.tra[0], s.tra[1], s.tra[2]); glTranslatef(s.tra[0], s.tra[1], s.tra[2]);
glScalef(s.sca, s.sca, s.sca); glScalef(s.sca, s.sca, s.sca);
float alpha; float alpha;
Point3f axis; Point3f axis;
s.rot.ToAxis(alpha, axis); s.rot.ToAxis(alpha, axis);
glRotatef(math::ToDeg(alpha), axis[0], axis[1], axis[2]); glRotatef(math::ToDeg(alpha), axis[0], axis[1], axis[2]);
} }
inline void glMultMatrix(const Similarityd &s) { inline void glMultMatrix(const Similarityd &s) {
@ -119,26 +120,27 @@ inline void glMultMatrix(const Similarityd &s) {
glScaled(s.sca, s.sca, s.sca); glScaled(s.sca, s.sca, s.sca);
} }
inline void glGetv(const GLenum pname, Matrix44f & m){ inline void glGetv(const GLenum pname, Eigen::Matrix<float,4,4>& matrix){
Matrix44f tmp; glGetFloatv(pname,matrix.data());
glGetFloatv(pname,tmp.V());
m = tmp.transpose();
} }
inline void glGetv(const GLenum pname, Eigen::Matrix<double,4,4>& matrix){
inline void glGetv(const GLenum pname, Matrix44d & m){ glGetDoublev(pname,matrix.data());
Matrix44d tmp;
glGetDoublev(pname,tmp.V());
m = tmp.transpose();
} }
inline void glGetv(const GLenum pname, Eigen::Matrix<float,4,4,Eigen::RowMajor>& matrix){
inline void glGetDirectv(const GLenum pname, Matrix44f & m){ Eigen::Matrix4f tmp;
glGetFloatv(pname,m.V()); glGetFloatv(pname,tmp.data());
matrix = tmp;
} }
inline void glGetv(const GLenum pname, Eigen::Matrix<double,4,4,Eigen::RowMajor>& matrix){
inline void glGetDirecv(const GLenum pname, Matrix44d & m){ Eigen::Matrix4d tmp;
glGetDoublev(pname,m.V()); glGetDoublev(pname,tmp.data());
matrix = tmp;
} }
template<typename Scalar>
inline void glGetv(const GLenum pname, const Eigen::Transform<Scalar,3>& t)
{ glGetv(pname,t.matrix()); }
}//namespace }//namespace
#endif #endif
#endif

68
wrap/gui/trackball.cpp
View File

@ -8,7 +8,7 @@
* \ * * \ *
* All rights reserved. * * All rights reserved. *
* * * *
* This program is free software; you can redistribute it and/or modify * * This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by * * it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or * * the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. * * (at your option) any later version. *
@ -124,12 +124,12 @@ Transform::Transform() {
} }
Trackball::Trackball(): current_button(0), current_mode(NULL), inactive_mode(NULL), Trackball::Trackball(): current_button(0), current_mode(NULL), inactive_mode(NULL),
dragging(false), spinnable(true), spinning(false), dragging(false), spinnable(true), spinning(false),
history_size(10){ history_size(10){
setDefaultMapping (); setDefaultMapping ();
} }
Trackball::~Trackball() Trackball::~Trackball()
{ {
std::map<int, TrackMode *>::iterator it; std::map<int, TrackMode *>::iterator it;
for(it = modes.begin(); it != modes.end(); it++) for(it = modes.begin(); it != modes.end(); it++)
@ -167,7 +167,7 @@ void Trackball::GetView() {
} }
// the drawing code has been moved to the trackmodes // the drawing code has been moved to the trackmodes
void Trackball::DrawPostApply() { void Trackball::DrawPostApply() {
if(current_mode !=NULL){ if(current_mode !=NULL){
current_mode->Draw(this); current_mode->Draw(this);
}else{ }else{
@ -177,12 +177,12 @@ void Trackball::DrawPostApply() {
} }
void Trackball::Apply () { void Trackball::Apply () {
glTranslate (center); glTranslate (center);track.Matrix().print();
glMultMatrix (track.Matrix()); glMultMatrix (track.Matrix());
glTranslate (-center); glTranslate (-center);
} }
void Trackball::Apply(bool ToDraw) { void Trackball::Apply(bool ToDraw) {
Apply(); Apply();
if(ToDraw){ if(ToDraw){
DrawPostApply(); DrawPostApply();
@ -190,7 +190,7 @@ void Trackball::Apply(bool ToDraw) {
} }
void Trackball::ApplyInverse() { void Trackball::ApplyInverse() {
glTranslate(center); glTranslate(center);
glMultMatrix(track.InverseMatrix()); glMultMatrix(track.InverseMatrix());
glTranslate(-center); glTranslate(-center);
@ -242,15 +242,15 @@ void Trackball::DrawPlane() {
void Trackball::ToAscii(char* result){ void Trackball::ToAscii(char* result){
float * f = (float*) &track; float * f = (float*) &track;
sprintf(result, "trackball(%f,%f,%f,%f,%f,%f,%f,%f,%f)", sprintf(result, "trackball(%f,%f,%f,%f,%f,%f,%f,%f,%f)",
f[0],f[1],f[2],f[3],f[4],f[5],f[6],f[7],f[8] ); f[0],f[1],f[2],f[3],f[4],f[5],f[6],f[7],f[8] );
} }
bool Trackball::SetFromAscii(const char * st){ bool Trackball::SetFromAscii(const char * st){
float * f = (float*) &track; float * f = (float*) &track;
int res= sscanf(st, "trackball(%f,%f,%f,%f,%f,%f,%f,%f,%f)", int res= sscanf(st, "trackball(%f,%f,%f,%f,%f,%f,%f,%f,%f)",
f+0,f+1,f+2,f+3,f+4,f+5,f+6,f+7,f+8 ); f+0,f+1,f+2,f+3,f+4,f+5,f+6,f+7,f+8 );
return (res==9); return (res==9);
} }
@ -279,7 +279,7 @@ void Trackball::DrawPlaneHandle() {
void Trackball::DrawIcon() { void Trackball::DrawIcon() {
glPushMatrix(); glPushMatrix();
glScale(radius); glScale(radius);
/// Here start the real drawing stuff /// Here start the real drawing stuff
float amb[4] ={.3f,.3f,.3f,1.0f}; float amb[4] ={.3f,.3f,.3f,1.0f};
@ -287,31 +287,31 @@ void Trackball::DrawIcon() {
//float col2[4]={.9f,.9f,1.0f,1.0f}; //float col2[4]={.9f,.9f,1.0f,1.0f};
glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_CURRENT_BIT | GL_LIGHTING_BIT); glPushAttrib(GL_ENABLE_BIT | GL_LINE_BIT | GL_CURRENT_BIT | GL_LIGHTING_BIT);
if(current_mode == NULL ) glLineWidth(DH.LineWidthStill); if(current_mode == NULL ) glLineWidth(DH.LineWidthStill);
else glLineWidth(DH.LineWidthMoving); else glLineWidth(DH.LineWidthMoving);
glEnable(GL_LIGHTING); glEnable(GL_LIGHTING);
glEnable(GL_LIGHT0); glEnable(GL_LIGHT0);
glEnable(GL_LINE_SMOOTH); glEnable(GL_LINE_SMOOTH);
glEnable(GL_BLEND); glEnable(GL_BLEND);
glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA); glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
glColor(DH.color); glColor(DH.color);
glMaterialfv(GL_FRONT_AND_BACK,GL_EMISSION,amb); glMaterialfv(GL_FRONT_AND_BACK,GL_EMISSION,amb);
glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,col); glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,col);
glPushMatrix(); glPushMatrix();
DrawCircle(); DrawCircle();
glPushMatrix(); glPushMatrix();
glRotatef(90,1,0,0); glRotatef(90,1,0,0);
DrawCircle(); DrawCircle();
glRotatef(90,0,1,0); glRotatef(90,0,1,0);
DrawCircle(); DrawCircle();
glPopMatrix(); glPopMatrix();
glPopMatrix(); glPopMatrix();
//glColor4f(1.0,.8f,.8f,1.0f); //glColor4f(1.0,.8f,.8f,1.0f);
glPopAttrib(); glPopAttrib();
@ -336,36 +336,36 @@ void Trackball::Reset() {
//interface //interface
void Trackball::MouseDown(int button) { void Trackball::MouseDown(int button) {
undo_track = track; undo_track = track;
current_button |= button; current_button |= button;
SetCurrentAction(); SetCurrentAction();
Hits.clear(); Hits.clear();
} }
void Trackball::MouseDown(int x, int y, int button) { void Trackball::MouseDown(int x, int y, int button) {
undo_track = track; undo_track = track;
current_button |= button; current_button |= button;
SetCurrentAction(); SetCurrentAction();
last_point = Point3f((float)x, (float)y, 0); last_point = Point3f((float)x, (float)y, 0);
Hits.clear(); Hits.clear();
} }
void Trackball::MouseMove(int x, int y) { void Trackball::MouseMove(int x, int y) {
if(current_mode == NULL) return; if(current_mode == NULL) return;
if(last_point[2] == -1) { //changed mode in the middle of moving if(last_point[2] == -1) { //changed mode in the middle of moving
last_point = Point3f((float)x, (float)y, 0); last_point = Point3f((float)x, (float)y, 0);
return; return;
} }
undo_track = track; undo_track = track;
current_mode->Apply(this, Point3f(float(x), float(y), 0)); current_mode->Apply(this, Point3f(float(x), float(y), 0));
} }
void Trackball::MouseUp(int /* x */, int /* y */, int button) { void Trackball::MouseUp(int /* x */, int /* y */, int button) {
undo_track = track; undo_track = track;
current_button &= (~button); current_button &= (~button);
SetCurrentAction(); SetCurrentAction();
} }
// it assumes that a notch of 1.0 is a single step of the wheel // it assumes that a notch of 1.0 is a single step of the wheel
void Trackball::MouseWheel(float notch) void Trackball::MouseWheel(float notch)
{ {
undo_track = track; undo_track = track;
int buttons = current_button; int buttons = current_button;
@ -375,7 +375,7 @@ void Trackball::MouseWheel(float notch)
{ {
ScaleMode scalemode; ScaleMode scalemode;
scalemode.Apply (this, notch); scalemode.Apply (this, notch);
} }
else else
{ {
current_mode->Apply(this, notch); current_mode->Apply(this, notch);
@ -405,7 +405,7 @@ void Trackball::ButtonDown(Trackball::Button button) {
if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) { if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) {
old_sticky = modes[current_button]->isSticky(); old_sticky = modes[current_button]->isSticky();
} }
current_button |= button; current_button |= button;
if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) { if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) {
new_sticky = modes[current_button]->isSticky(); new_sticky = modes[current_button]->isSticky();
} }
@ -414,13 +414,13 @@ void Trackball::ButtonDown(Trackball::Button button) {
SetCurrentAction(); SetCurrentAction();
} }
void Trackball::ButtonUp(Trackball::Button button) { void Trackball::ButtonUp(Trackball::Button button) {
bool old_sticky=false, new_sticky=false; bool old_sticky=false, new_sticky=false;
assert ( modes.count (0) ); assert ( modes.count (0) );
if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) { if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) {
old_sticky = modes[current_button]->isSticky(); old_sticky = modes[current_button]->isSticky();
} }
current_button &= (~button); current_button &= (~button);
if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) { if ( ( modes.count (current_button) ) && ( modes[current_button] != NULL ) ) {
new_sticky = modes[current_button]->isSticky(); new_sticky = modes[current_button]->isSticky();
} }
@ -440,12 +440,12 @@ void Trackball::Undo(){
void Trackball::SetSpinnable(bool /* on*/ ){} void Trackball::SetSpinnable(bool /* on*/ ){}
bool Trackball::IsSpinnable() { bool Trackball::IsSpinnable() {
return spinnable; return spinnable;
} }
void Trackball::SetSpinning(Quaternionf &/* spin*/){} void Trackball::SetSpinning(Quaternionf &/* spin*/){}
void Trackball::StopSpinning(){} void Trackball::StopSpinning(){}
bool Trackball::IsSpinning() { bool Trackball::IsSpinning() {
return spinning; return spinning;
} }
//navigation interface: //navigation interface:
void Trackball::Back(){} void Trackball::Back(){}
@ -469,14 +469,14 @@ void Trackball::SetCurrentAction ()
} }
////return center of trackball in Window coordinates. ////return center of trackball in Window coordinates.
//Point3f Trackball::ScreenOrigin() { //Point3f Trackball::ScreenOrigin() {
// return camera.Project(ModelOrigin()); // return camera.Project(ModelOrigin());
//} //}
//return center of trackball in Model coordinates //return center of trackball in Model coordinates
//Point3f Trackball::ModelOrigin() { //Point3f Trackball::ModelOrigin() {
// return center; // return center;
//} //}
//Matrix44f Trackball::ScreenToModel() { //Matrix44f Trackball::ScreenToModel() {