Click here to Skip to main content
15,894,720 members
Articles / Desktop Programming / WTL

Part One: Euclidian Geometric Algebra and Quaternions

Rate me:
Please Sign up or sign in to vote.
4.87/5 (12 votes)
15 Mar 2005CPOL8 min read 94.7K   2.3K   63  
Geometric Algebra applied to OpenGL
/*
Gaigen, Copyright (c) 2001-2004, University of Amsterdam

Copying, use and development for education and research purposes
permitted as long as this license is not removed from the files.

All rights for commercial use reserved; for more information
contact Daniel Fontijne (fontijne@science.uva.nl)

This software is unsupported.
*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#ifdef GAIM_NAMESPACE
namespace CLASSNAME_NS {
#endif // GAIM_NAMESPACE

CLASSNAME::CLASSNAME() {

}
/*
CLASSNAME::CLASSNAME(GAIM_CLASSNAME &a) {
	copy(a);
}
*/
CLASSNAME::~CLASSNAME() {

}

GAIM_RETURN_TYPE CLASSNAME::operator=(const CLASSNAME &a) {
	copy(a);
	return *this;
}

GAIM_RETURN_TYPE CLASSNAME::operator=(GAIM_FLOAT f) {
	set(GRADE0, &f);
	return *this;
}

#ifdef GAIM_FUNCTION_FASTTEMPVAR
GAIM_RETURN_TYPE CLASSNAME::getTemp() {
	// gets a temporary variable (based on Stroustup (heaven help you if you write an expression which requires more than MV_MAX_TEMP temporaries))
	static int idx = 0;
	static CLASSNAME *temp = NULL;
	if (!temp) {
		int i;
		temp = new CLASSNAME[MV_MAX_TEMP];
		for (i = 0; i < MV_MAX_TEMP; i++) {
			temp[i].setUsage(0xf | (0x7 << 4));
			temp[i].setUsage(0);
		}
	}
//	printf("Temp idx = %d\n", idx);
	if (idx == MV_MAX_TEMP) idx = 0;
	return temp[idx++];
}
#endif // GAIM_FUNCTION_FASTTEMPVAR

#ifdef GAIM_FUNCTION_TAKEGRADE
GAIM_RETURN_TYPE CLASSNAME::grade(int g) const {
	GAIM_RETURN_VAR(result);
//	GAIM_RETURN_VAR(result);
	result.takeGrade(*this, g);
	return result;
}
#endif // GAIM_FUNCTION_TAKEGRADE


#ifdef GAIM_FUNCTION_REVERSE
GAIM_RETURN_TYPE CLASSNAME::operator~() const {
	GAIM_RETURN_VAR(result);
	((GAIM_CLASSNAME&)result).reverse(*this);
	return result;
}
#endif // GAIM_FUNCTION_REVERSE

#ifdef GAIM_FUNCTION_CLIFFORDCONJUGATE
GAIM_RETURN_TYPE CLASSNAME::operator--() const {
	GAIM_RETURN_VAR(result);
	((GAIM_CLASSNAME&)result).cliffordConjugate(*this);
	return result;
}
#endif // GAIM_FUNCTION_CLIFFORDCONJUGATE

#ifdef GAIM_FUNCTION_GRADEINVOLUTION
GAIM_RETURN_TYPE CLASSNAME::operator++() const {
	GAIM_RETURN_VAR(result);
	((GAIM_CLASSNAME&)result).gradeInvolution(*this);
	return result;
}
#endif // GAIM_FUNCTION_GRADEINVOLUTION

#ifdef GAIM_PRODUCT_GP
GAIM_RETURN_TYPE CLASSNAME::operator*=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.gp(*this, a);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator*=(GAIM_FLOAT a) {
	CLASSNAME tmp;
	tmp.op(a, *this);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC gp(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.gp(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC gp(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.op(a, b);
	return result;
}

#define USE_NEW_EXP

#ifdef USE_OLD_EXP
GAIM_RETURN_TYPE CLASSNAME::exp(int order /*= 9*/) const {
	int i, div = 1;
	GAIM_RETURN_VAR(result);

	result.setScalar(1.0);
	if (order == 0) return result;

	CLASSNAME tmp(*this);

	for (i = 1; i < order; i++) {
		result += tmp * (GAIM_FLOAT)(1.0 / (GAIM_FLOAT)div);

		if (i < order-1) {
			div *= (i+1);
			tmp *= *this;
		}
	}

	return result;
}
#endif /* USE_OLD_EXP */

#ifdef USE_NEW_EXP
GAIM_RETURN_TYPE CLASSNAME::exp(int order /*= 9*/) const {
	/*
	Improved version of exp thanks to Robert Valkenburg & students
	*/
    int i;
    CLASSNAME result;

    result.setScalar(1.0);
    if (order == 0) {
        GAIM_RETURN_VAR(r);
        r = result;
        return r;
    }

    // scale by power of 2 so that its norm is < 1
    unsigned long max = (unsigned long)largestCoordinate();
    unsigned long scale=1;
    if (max > 1) scale <<= 1;
    while (max)
    {
        max >>= 1;
        scale <<= 1;
    }

    CLASSNAME scaled = (*this) / (GAIM_FLOAT)(scale);

    // taylor approximation
    CLASSNAME tmp;

    tmp.setScalar(1.0);
    for (i = 1; i < order; i++) {
        tmp = tmp*scaled/(GAIM_FLOAT)(i);
        result += tmp;
    }

    // undo scaling
    while (scale > 1)
    {
        result *= result;
        scale >>= 1;
    }
    GAIM_RETURN_VAR(r);
    r = result;
    return r;
}
#endif /* USE_NEW_EXP */


#ifndef GAIM_FUNCTION_FASTDUAL
GAIM_RETURN_TYPE CLASSNAME::dual() const {
	GAIM_RETURN_VAR(result);
	result.gp(*this, Ii);
	return result;
}
#endif // GAIM_FUNCTION_FASTDUAL

#endif // GAIM_PRODUCT_GP

#ifdef GAIM_PRODUCT_GP_EM
GAIM_RETURN_TYPE CLASSNAME_NS_DSC gpem (const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result. GAIM_PRODUCT_GP_EM (a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC gpem (GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.op(a, b);
	return result;
}
#endif // GAIM_PRODUCT_GP_EM


#ifdef GAIM_FUNCTION_FASTDUAL
GAIM_RETURN_TYPE CLASSNAME::fastDual() const {
	GAIM_RETURN_VAR(result);
	((GAIM_CLASSNAME&)result).fastDual(*this);
	return result;
}
#endif // GAIM_FUNCTION_FASTDUAL


#ifdef GAIM_FUNCTION_VERSORINVERSE
GAIM_RETURN_TYPE CLASSNAME::versorInverse() const {
	GAIM_RETURN_VAR(result);
	result.versorInverse(*this);
	return result;
}
#endif
#ifdef GAIM_FUNCTION_LOUNESTOINVERSE
GAIM_RETURN_TYPE CLASSNAME::lounestoInverse() const {
	GAIM_RETURN_VAR(result);
	result.lounestoInverse(*this);
	return result;
}
#endif
#ifdef GAIM_FUNCTION_GENERALINVERSE
GAIM_RETURN_TYPE CLASSNAME::generalInverse() const {
	GAIM_RETURN_VAR(result);
	result.generalInverse(*this);
	return result;
}
#endif

#ifdef GAIM_PRODUCT_IGP
GAIM_RETURN_TYPE CLASSNAME::inverse() const {
	GAIM_RETURN_VAR(result);
	result.inverse(*this);
	return result;
}

GAIM_RETURN_TYPE CLASSNAME::operator/=(const CLASSNAME &a) {
	CLASSNAME tmp, tmp2;
	tmp2.inverse(a);
	tmp.gp(*this, tmp2);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator/=(GAIM_FLOAT a) {
	CLASSNAME tmp;
	tmp.op(1.0f / a, *this);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_PRODUCT_IGP

#ifdef GAIM_PRODUCT_SCP
GAIM_RETURN_TYPE CLASSNAME::operator%=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.scp(*this, a);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator%=(GAIM_FLOAT a) {
	GAIM_FLOAT f = scalar() * a;
	setScalar(&f);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC scp(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.scp(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC scp(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.set(GRADE0, mulScalar(a, b.scalar()));
	return result;
}
#endif // GAIM_PRODUCT_SCP

#ifdef GAIM_PRODUCT_LCONT
GAIM_RETURN_TYPE CLASSNAME_NS_DSC lcont(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.lcont(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME::operator<<=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.lcont(*this, a);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_PRODUCT_LCONT

#ifdef GAIM_PRODUCT_LCONT_EM
GAIM_RETURN_TYPE CLASSNAME_NS_DSC lcem(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result. GAIM_PRODUCT_LCONT_EM (a, b);
	return result;
}
#endif // GAIM_PRODUCT_LCONT_EM

#ifdef GAIM_PRODUCT_RCONT
GAIM_RETURN_TYPE CLASSNAME_NS_DSC rcont(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.rcont(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME::operator>>=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.rcont(*this, a);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_PRODUCT_RCONT

#ifdef GAIM_PRODUCT_HIP
GAIM_RETURN_TYPE CLASSNAME_NS_DSC hip(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.hip(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC hip(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	CLASSNAME _a(a);
	result.hip(_a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC hip(const CLASSNAME &a, GAIM_FLOAT b) {
	GAIM_RETURN_VAR(result);
	CLASSNAME _b(b);
	result.hip(a, _b);
	return result;
}
#endif // GAIM_PRODUCT_HIP

#ifdef GAIM_PRODUCT_MHIP
GAIM_RETURN_TYPE CLASSNAME_NS_DSC mhip(const CLASSNAME &a, const CLASSNAME &b) { // GAIM_RETURN_TYPE used to read CLASSNAME & (modified 10-10-2002)
	GAIM_RETURN_VAR(result);
	result.mhip(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC mhip(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	CLASSNAME _a(a);
	result.mhip(_a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC mhip(const CLASSNAME &a, GAIM_FLOAT b) {
	GAIM_RETURN_VAR(result);
	CLASSNAME _b(b);
	result.mhip(a, _b);
	return result;
}
#endif // GAIM_PRODUCT_MHIP


#ifdef GAIM_PRODUCT_OP
GAIM_RETURN_TYPE CLASSNAME_NS_DSC op(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.op(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC op(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.op(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME::operator^=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.op(*this, a);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator^=(GAIM_FLOAT a) {
	CLASSNAME tmp;
	tmp.op(a, *this);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_PRODUCT_OP

#ifdef GAIM_FUNCTION_MEETJOIN
GAIM_RETURN_TYPE CLASSNAME_NS_DSC meet(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.meet(a, b);
	return result;
}

GAIM_RETURN_TYPE CLASSNAME_NS_DSC join(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.join(a, b);
	return result;
}

#endif // GAIM_FUNCTION_MEETJOIN


#ifdef GAIM_FUNCTION_ADD
GAIM_RETURN_TYPE CLASSNAME_NS_DSC add(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.add(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC add(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.add(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME::operator+=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.add(*this, a);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator+=(GAIM_FLOAT a) {
	CLASSNAME tmp;
	tmp.add(a, *this);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_FUNCTION_ADD

#ifdef GAIM_FUNCTION_SUBSTRACT
GAIM_RETURN_TYPE CLASSNAME_NS_DSC sub(const CLASSNAME &a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.sub(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME_NS_DSC sub(GAIM_FLOAT a, const CLASSNAME &b) {
	GAIM_RETURN_VAR(result);
	result.sub(a, b);
	return result;
}
GAIM_RETURN_TYPE CLASSNAME::operator-=(const CLASSNAME &a) {
	CLASSNAME tmp;
	tmp.sub(*this, a);
	this->copy(tmp);
	return *this;
}
GAIM_RETURN_TYPE CLASSNAME::operator-=(GAIM_FLOAT a) {
	CLASSNAME tmp;
	tmp.add(a, *this);
	this->copy(tmp);
	return *this;
}
#endif // GAIM_FUNCTION_SUBSTRACT

#ifdef GAIM_FUNCTION_NEGATE
GAIM_RETURN_TYPE CLASSNAME_NS_DSC negate(const CLASSNAME &a) {
	GAIM_RETURN_VAR(result);
	((GAIM_CLASSNAME&)result).negate(a);
	return result;
}
#endif //GAIM_FUNCTION_NEGATE



#ifdef GAIM_FUNCTION_NORMALIZE
GAIM_RETURN_TYPE CLASSNAME::normal(int norm /* = 1 */) const {
	GAIM_RETURN_VAR(result);
	result.normalize(*this, norm);
	return result;
}
#endif // GAIM_FUNCTION_NORMALIZE


#ifdef GAIM_PRODUCT_OM
CLASSNAME_OM::CLASSNAME_OM(const CLASSNAME *vectorImages[3]) {
	initVectorImages(vectorImages);
}

CLASSNAME_OM::CLASSNAME_OM(const CLASSNAME &spinor) {
	initSpinor(spinor);
}


int CLASSNAME_OM::initSpinor(const CLASSNAME &spinor) {
	CLASSNAME si(spinor.inverse()), v, vi[GA_MAX_DIM];
	GAIM_FLOAT coordinates[GA_MAX_DIM];
	int i;

	memset(coordinates, 0, sizeof(GAIM_FLOAT) * GA_MAX_DIM);

	for (i = 0; i < CLASSNAME::dim; i++) {
		coordinates[i] = 1.0;
		if (i > 0) coordinates[i-1] = 0.0;
		v.setVector(coordinates);
		vi[i] = ((spinor * v) * si)(GRADE1);

	}
	return initVectorImages(vi);
}


GAIM_RETURN_TYPE CLASSNAME_NS_DSC om(const CLASSNAME_OM &om, const CLASSNAME &a) {
	GAIM_RETURN_VAR(result);

	result.om(a, om);

	return result;
}


#ifdef GAIM_NAMESPACE
}
#endif // GAIM_NAMESPACE


#endif // GAIM_PRODUCT_OM


By viewing downloads associated with this article you agree to the Terms of Service and the article's licence.

If a file you wish to view isn't highlighted, and is a text file (not binary), please let us know and we'll add colourisation support for it.

License

This article, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)


Written By
Systems Engineer
United States United States
In 2001 my wife and I started on a new business venture. Reserve Analyst Software

http://www.ReserveAnalyst.com

I have been messing with computers since the 6502 was announced. A good deal on the hardware side.


http://www.Lakeweb.net


Thanks, Dan.


Comments and Discussions