dvec.h
Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #ifndef __DVEC_H__
00026 #define __DVEC_H__
00027
00028 #include "common.h"
00029
00030 #include <math.h>
00031
00032 #include "raytrace.h"
00033
00034
00035 extern "C++" {
00036 #include <iostream>
00037
00038 const double VEQUALITY = 0.0000001;
00039
00040 template<int LEN>
00041 struct vec_internal;
00042
00043 template<int LEN>
00044 class dvec;
00045
00046 template<int LEN>
00047 std::ostream& operator<<(std::ostream& out, const dvec<LEN>& v);
00048
00049 class dvec_unop {
00050 public:
00051 virtual double operator()(double a) const = 0;
00052 virtual ~dvec_unop() {}
00053 };
00054
00055 class dvec_op {
00056 public:
00057 virtual double operator()(double a, double b) const = 0;
00058 virtual ~dvec_op() {}
00059 };
00060
00061 template<int LEN>
00062 class dvec {
00063 public:
00064 dvec(double s);
00065 dvec(const double* vals, bool aligned=true);
00066 dvec(const dvec<LEN>& p);
00067
00068 dvec<LEN>& operator=(const dvec<LEN>& p);
00069 double operator[](int index) const;
00070 void u_store(double* arr) const;
00071 void a_store(double* arr) const;
00072
00073 bool operator==(const dvec<LEN>& b) const;
00074
00075 dvec<LEN> operator+(const dvec<LEN>& b);
00076 dvec<LEN> operator-(const dvec<LEN>& b);
00077 dvec<LEN> operator*(const dvec<LEN>& b);
00078 dvec<LEN> operator/(const dvec<LEN>& b);
00079
00080 dvec<LEN> madd(const dvec<LEN>& s, const dvec<LEN>& b);
00081 dvec<LEN> madd(const double s, const dvec<LEN>& b);
00082
00083 dvec<LEN> map(const dvec_unop& operation, int limit = LEN);
00084 double foldr(double proto, const dvec_op& operation, int limit = LEN);
00085 double foldl(double proto, const dvec_op& operation, int limit = LEN);
00086
00087 class mul : public dvec_op {
00088 public:
00089 double operator()(double a, double b) const { return a * b; }
00090 };
00091
00092 class add : public dvec_op {
00093 public:
00094 double operator()(double a, double b) const { return a + b; }
00095 };
00096
00097 class sub : public dvec_op {
00098 public:
00099 double operator()(double a, double b) const { return a - b; }
00100 };
00101
00102 class sqrt : public dvec_unop {
00103 public:
00104 double operator()(double a) const { return ::sqrt(a); }
00105 };
00106 private:
00107 vec_internal<LEN> data;
00108
00109 dvec(const vec_internal<LEN>& d);
00110 };
00111
00112
00113
00114
00115 #define VEC_ALIGN
00116
00117
00118 #if defined(__SSE2__) && defined(__GNUC__) && defined(HAVE_EMMINTRIN_H)
00119 # define __x86_vector__
00120 # include "vector_x86.h"
00121 #else
00122 # define __fpu_vector__
00123 # include "vector_fpu.h"
00124 #endif
00125
00126 inline bool vequals(const vec2d& a, const vec2d& b) {
00127 return
00128 (fabs(a.x()-b.x()) < VEQUALITY) &&
00129 (fabs(a.y()-b.y()) < VEQUALITY);
00130 }
00131
00132
00133 typedef fastf_t mat2d_t[4] VEC_ALIGN;
00134
00135
00136 typedef fastf_t pt2d_t[2] VEC_ALIGN;
00137
00138
00139
00140 inline
00141 bool mat2d_inverse(mat2d_t inv, mat2d_t m) {
00142 pt2d_t _a = {m[0], m[1]};
00143 pt2d_t _b = {m[3], m[2]};
00144 dvec<2> a(_a);
00145 dvec<2> b(_b);
00146 dvec<2> c = a*b;
00147 fastf_t det = c.foldr(0, dvec<2>::sub());
00148 if (NEAR_ZERO(det, VUNITIZE_TOL)) return false;
00149 fastf_t scale = 1.0 / det;
00150 double tmp[4] VEC_ALIGN = {m[3], -m[1], -m[2], m[0]};
00151 dvec<4> iv(tmp);
00152 dvec<4> sv(scale);
00153 dvec<4> r = iv * sv;
00154 r.a_store(inv);
00155 return true;
00156 }
00157 inline
00158 void mat2d_pt2d_mul(pt2d_t r, mat2d_t m, pt2d_t p) {
00159 pt2d_t _a = {m[0], m[2]};
00160 pt2d_t _b = {m[1], m[3]};
00161 dvec<2> x(p[0]);
00162 dvec<2> y(p[1]);
00163 dvec<2> a(_a);
00164 dvec<2> b(_b);
00165 dvec<2> c = a*x + b*y;
00166 c.a_store(r);
00167 }
00168 inline
00169 void pt2dsub(pt2d_t r, pt2d_t a, pt2d_t b) {
00170 dvec<2> va(a);
00171 dvec<2> vb(b);
00172 dvec<2> vr = va - vb;
00173 vr.a_store(r);
00174 }
00175
00176 inline
00177 fastf_t v2mag(pt2d_t p) {
00178 dvec<2> a(p);
00179 dvec<2> sq = a*a;
00180 return sqrt(sq.foldr(0, dvec<2>::add()));
00181 }
00182 inline
00183 void move(pt2d_t a, const pt2d_t b) {
00184 a[0] = b[0];
00185 a[1] = b[1];
00186 }
00187 }
00188
00189 #endif
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199