home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
DOS/V Power Report 1997 March
/
VPR9703A.ISO
/
VPR_DATA
/
DOGA
/
SOURCES
/
PASM.LZH
/
MATRIX.CPP
< prev
next >
Wrap
C/C++ Source or Header
|
1996-05-27
|
6KB
|
278 lines
#include <math.h>
#include "matrix.h"
#include "log.h"
Matrix::Matrix(int n)
{
v[0] = v[1] = v[2] = v[3] = Vector(0.0,0.0,0.0);
if (n != 0) {
v[0].x = v[1].y = v[2].z = 1.0;
}
}
Matrix::Matrix(const Matrix& m)
{
v[0]=m.v[0];
v[1]=m.v[1];
v[2]=m.v[2];
v[3]=m.v[3];
}
Matrix::Matrix(const Vector& v0, const Vector& v1, const Vector& v2, const Vector& v3)
{
v[0] = v0;
v[1] = v1;
v[2] = v2;
v[3] = v3;
}
Matrix::Matrix(const Vector& vx, const Vector& vz)
{
*this = vx ^ vz;
}
Matrix& Matrix::operator +=(const Matrix& m)
{
v[0] += m.v[0];
v[1] += m.v[1];
v[2] += m.v[2];
v[3] += m.v[3];
return *this;
}
Matrix& Matrix::operator -=(const Matrix& m)
{
v[0] -= m.v[0];
v[1] -= m.v[1];
v[2] -= m.v[2];
v[3] -= m.v[3];
return *this;
}
const int Matrix::operator ==(const Matrix& m2) const
{
return v[0] == m2.v[0]
&& v[1] == m2.v[1]
&& v[2] == m2.v[2]
&& v[3] == m2.v[3];
}
const int Matrix::operator !=(const Matrix& m2) const
{
return !(*this == m2);
}
const Matrix Matrix::operator +(const Matrix& m2) const
{
Matrix temp(*this);
temp += m2;
return temp;
}
const Matrix Matrix::operator -(const Matrix& m2) const
{
Matrix temp(*this);
temp -= m2;
return temp;
}
const Vector Matrix::operator *(const Vector& v1) const
{
Vector temp;
temp[0] = v[0].x * v1.x + v[1].x * v1.y + v[2].x * v1.z + v[3].x;
temp[1] = v[0].y * v1.x + v[1].y * v1.y + v[2].y * v1.z + v[3].y;
temp[2] = v[0].z * v1.x + v[1].z * v1.y + v[2].z * v1.z + v[3].z;
return temp;
}
const Vector Matrix::mul_notmove(const Vector &v1) const
{
Vector temp;
temp.x = v[0].x * v1.x + v[1].x * v1.y + v[2].x * v1.z;
temp.y = v[0].y * v1.x + v[1].y * v1.y + v[2].y * v1.z;
temp.z = v[0].z * v1.x + v[1].z * v1.y + v[2].z * v1.z;
return temp;
}
const Matrix Matrix::operator *(const Matrix& m2) const
{
Matrix temp;
temp.v[0] = mul_notmove(m2.v[0]);
temp.v[1] = mul_notmove(m2.v[1]);
temp.v[2] = mul_notmove(m2.v[2]);
temp.v[3] = mul_notmove(m2.v[3]) + v[3];
return temp;
}
const Matrix Matrix::m_move(const Vector& v1)
{
Matrix temp(1);
temp.v[3] = v1;
return temp;
}
const Matrix Matrix::m_scale(const Vector& v1)
{
Matrix temp(0);
temp.v[0].x = v1.x;
temp.v[1].y = v1.y;
temp.v[2].z = v1.z;
return temp;
}
const Matrix Matrix::m_rotx(double t)
{
Matrix temp(1);
temp.v[1].y = cos(t);
temp.v[1].z = sin(t);
temp.v[2].y = -temp.v[1][2];
temp.v[2].z = temp.v[1][1];
return temp;
}
const Matrix Matrix::m_roty(double t)
{
Matrix temp(1);
temp.v[2].z = cos(t);
temp.v[2].x = sin(t);
temp.v[0].z = -temp.v[2].x;
temp.v[0].x = temp.v[2].z;
return temp;
}
const Matrix Matrix::m_rotz(double t)
{
Matrix temp(1);
temp.v[0].x = cos(t);
temp.v[0].y = sin(t);
temp.v[1].x = -temp.v[0].y;
temp.v[1].y = temp.v[0].x;
return temp;
}
const double Matrix::value(void) const
{
return
v[0].x*v[1].y*v[2].z+v[1].x*v[2].y*v[0].z+v[2].x*v[0].y*v[1].z
-v[0].x*v[2].y*v[1].z-v[1].x*v[0].y*v[2].z-v[2].x*v[1].y*v[0].z;
}
static double value(const Vector& v0, const Vector& v1, const Vector& v2)
{
return v0.x*v1.y*v2.z+v1.x*v2.y*v0.z+v2.x*v0.y*v1.z
-v0.x*v2.y*v1.z-v1.x*v0.y*v2.z-v2.x*v1.y*v0.z;
}
const Matrix Matrix::inv(void) const
{
Matrix dst(0);
double val = value();
if (-minimumdouble < val && val < minimumdouble) {
return Matrix(1);
}
dst.v[0].x = (v[1].y*v[2].z-v[2].y*v[1].z)/val;
dst.v[0].y = -(v[0].y*v[2].z-v[2].y*v[0].z)/val;
dst.v[0].z = (v[0].y*v[1].z-v[1].y*v[0].z)/val;
dst.v[1].x = -(v[1].x*v[2].z-v[2].x*v[1].z)/val;
dst.v[1].y = (v[0].x*v[2].z-v[2].x*v[0].z)/val;
dst.v[1].z = -(v[0].x*v[1].z-v[1].x*v[0].z)/val;
dst.v[2].x = (v[1].x*v[2].y-v[2].x*v[1].y)/val;
dst.v[2].y = -(v[0].x*v[2].y-v[2].x*v[0].y)/val;
dst.v[2].z = (v[0].x*v[1].y-v[1].x*v[0].y)/val;
if (v[3].x != 0.0 || v[3].y != 0.0 || v[3].z != 0.0) {
dst.v[3].x = -::value(v[1], v[2], v[3])/val;
dst.v[3].y = ::value(v[0], v[2], v[3])/val;
dst.v[3].z = -::value(v[0], v[1], v[3])/val;
}
return dst;
}
const Matrix Matrix::tra(void) const
{
Matrix dst(0);
dst.v[0].x = v[0].x;
dst.v[0].y = v[1].x;
dst.v[0].z = v[2].x;
dst.v[1].x = v[0].y;
dst.v[1].y = v[1].y;
dst.v[1].z = v[2].y;
dst.v[2].x = v[0].z;
dst.v[2].y = v[1].z;
dst.v[2].z = v[2].z;
return dst;
}
const Vector Matrix::GetRotation(void) const
{
Vector rot;
if (v[0].z <= -1.0) {
rot.y = deg(90);
} else if (v[0].z >= 1.0) {
rot.y = deg(-90);
} else {
rot.y = asin(-v[0].z);
}
double s = 1.0 - v[0].z * v[0].z;
if (s <= minimumdouble) {
rot.x = 0.0;
if (-minimumdouble < v[1].x && v[1].x < minimumdouble
&& -minimumdouble < v[1].y && v[1].y < minimumdouble) {
rot.z = 0.0;
} else {
rot.z = atan2(-v[1].x, v[1].y);
}
} else {
double cosy = sqrt(s);
if (-minimumdouble < v[1].z && v[1].z < minimumdouble
&& -minimumdouble < v[2].z && v[2].z < minimumdouble) {
rot.x = 0.0;
} else {
rot.x = atan2(v[1].z/cosy, v[2].z/cosy);
}
if (-minimumdouble < v[0].y && v[0].y < minimumdouble
&& -minimumdouble < v[0].x && v[0].x < minimumdouble) {
rot.z = 0.0;
} else {
rot.z = atan2(v[0].y/cosy, v[0].x/cosy);
}
}
return rot;
}
const Matrix Vector::operator ^(const Vector& vz) const
{
Matrix m(0);
double l = length();
if (l <= minimumdouble) {
return Matrix(1);
}
m.v[0] = (*this) * (1.0 / l);
Vector vy = vz * (*this);
double vyl = vy.length();
if (vyl < minimumdouble) {
vy = Vector(0.0, 0.0, 1.0) * (*this);
vyl= vy.length();
if (vyl < minimumdouble) {
vy = Vector(1.0, 0.0, 0.0) * (*this);
vyl= vy.length();
if (vyl < minimumdouble) {
return Matrix(1);
}
}
}
m.v[1] = vy * (1.0 / vyl);
// Vector vz2 = (*this) * vy;
// m.v[2] = (1.0 / vz2.length()) * vz2;
m.v[2] = m.v[0] * m.v[1];
return m;
}