geometry.h

#ifndef GEOMETRY_H
#define GEOMETRY_H
/*
 * 3-D geometry (matrix, vector, quaternion) functions for partiview.
 *
 * Stuart Levy, slevy@ncsa.uiuc.edu
 * National Center for Supercomputing Applications,
 * University of Illinois 2001.
 * This file is part of partiview, released under the
 * Illinois Open Source License; see the file LICENSE.partiview for details.
 */

#ifdef __cplusplus
extern "C" {
#endif

#ifndef CONST
# ifdef __cplusplus
#  define CONST const
# else
#  define CONST
# endif
#endif

#ifdef unix
# include <alloca.h>
#endif

typedef struct { float x[3]; } Point;
typedef struct { float q[4]; } Quat;  /* {r,i,j,k} => {0,1,2,3} */

typedef struct { float m[4*4]; } Matrix;

#define NewA(type, count)  (type *)alloca((count) * sizeof(type))

#ifndef COUNT
# define COUNT(array)  (sizeof(array) / sizeof((array)[0]))
#endif

/*
 * Global Constants
 */

extern Matrix Tidentity;


/*
 * Global variables
 */

void arena_init( int nbytes );

void vsub( Point *dst, CONST Point *a, CONST Point *b );
void vadd( Point *dst, CONST Point *a, CONST Point *b );
void vcross( Point *dst, CONST Point *a, CONST Point *b );
float vdot( CONST Point *a, CONST Point *b );
void vscale( Point *dst, float s, CONST Point *src );
void vsadd( Point *dst, CONST Point *a, float sb, CONST Point *b );
void vlerp( Point *dst, float frac, CONST Point *vfrom, CONST Point *vto );
void vcomb( Point *dst, float sa, CONST Point *a, float sb, CONST Point *b );
float vunit( Point *dst, CONST Point *src );
void vproj( Point *along, Point *perp, CONST Point *vec, CONST Point *onto );

float qdot( CONST Quat *q1, CONST Quat *q2 );
void qnorm( Quat *qdst, CONST Quat *qsrc );
void qcomb( Quat *qdst, float sa, CONST Quat *qa, float sb, CONST Quat *qb );

        /* Distance between two 3-D points */
float vdist( CONST Point *p1, CONST Point *p2 );
        /* Distance between two quaternions */
float iqdist( CONST Point *q1, CONST Point *q2 );
float qdist( CONST Quat *q1, CONST Quat *q2 );
        /* Distance between two matrices */
float tdist( CONST Matrix *t1, CONST Matrix *t2 );

        /* vector magnitude */
float vlength( CONST Point *v );

        /* Transform a 3-D point by a matrix: dst = src*T. */
void vtfmpoint( Point *dst,  CONST Point *src, CONST Matrix *T );

        /* Transform a vector by a (Euclidean) matrix: dst = src*T */
void vtfmvector( Point *dst,  CONST Point *src, CONST Matrix *T );
void vuntfmvector( Point *dst,  CONST Point *src, CONST Matrix *T );

        /* Get translation part of matrix */
void vgettranslation( Point *dst, CONST Matrix *T );
        /* Set translation part of matrix, leave 3x3 submatrix alone */
void vsettranslation( Matrix *T, CONST Point *src );

void vrotxy( Point *dst, CONST Point *src, CONST float cs[2] );
void mcopy( Matrix *dst, CONST Matrix *src );

        /* 4x4 matrix multiply */
void mmmul( Matrix *dst, CONST Matrix *a, CONST Matrix *b );

        /* Matrix inverse for Euclidean similarities: matrices are
         * a product of pure translation/rotation/uniform-scaling --
         */
void eucinv( Matrix *dst, CONST Matrix *src );

        /* matrix conjugation: see its use in Gview.C */
void mconjugate( Matrix *To2wout, CONST Matrix *To2win, CONST Matrix *Tincrf,
                 CONST Matrix *Tf2w, CONST Matrix *Tw2f,
                 CONST Point *pcenw, CONST Point *pcenf );

void grotation( Matrix *Trot, CONST Point *fromaxis, CONST Point *toaxis );

        /* matrix to 3-component quaternion (real-part is sqrt(1-iquat.iquat),
         * and real part (in range 0 .. +1) is also returned as func value.
         */
float tfm2iquat( Point *iquat, CONST Matrix *src );     /* returns real part */
void tfm2quat( Quat *quat, CONST Matrix *src );
void iquat2tfm( Matrix *dst, CONST Point *iquat );
void quat2tfm( Matrix *dst, CONST Quat *quat );

        /* quaternion interpolation */
void iquat_lerp( Point *dquat, float frac, CONST Point *qfrom, CONST Point *qto );
void quat_lerp( Quat *dquat, float frac, CONST Quat *qfrom, CONST Quat *qto );

        /* axis/angle rotation -> quaternion */
void rot2quat( Quat *quat, float degrees, CONST Point *axis );
void rot2iquat( Point *iquat, float degrees, CONST Point *axis );

        /* construct matrix from axis/angle rotation */
void rot2tfm( Matrix *dst, float degrees, CONST Point *axis );

        /* Construct a matrix that's a pure X/Y/Z rotation/scaling/translation */
void mrotation( Matrix *dst, float degrees, char xyzaxis );
void mscaling( Matrix *dst, float sx, float sy, float sz );
void mtranslation( Matrix *dst, float tx, float ty, float tz );

        /* cam2world <-> az(Y) el(X) roll(Z) + translation */
float tfm2xyzaer( Point *xyz, float aer[3], CONST Matrix *c2w ); /* => scale */
void xyzaer2tfm( Matrix *c2w, CONST Point *xyz, CONST float aer[3] );

        /* perspective frustum */
void mfrustum( Matrix *Tproj, float nl, float nr, float nd, float nu, float cnear, float cfar );


#ifdef __cplusplus
}
#endif


#endif