OSDN Git Service

fix for build on winXP
[moflib/moflib.git] / src / mof / math / basic_matrix.hpp
index 6913abf..4213283 100644 (file)
@@ -2,6 +2,7 @@
 #include <mof/base/mofdef.hpp>
 #include <mof/math/threshold.hpp>
 #include <boost/operators.hpp>
+#include <boost/utility.hpp>
 #include <ostream>
 #include <iomanip>
 #include <cmath>
@@ -10,8 +11,11 @@ namespace mof
 {
 namespace math
 {
+       template <size_t Dim> class row_of_matrix;
+
        /**
-        * @brief アフィン同次変換行列テンプレートクラス
+        * @brief 同次座標変換行列テンプレートクラス
+        * @note  このクラスはアフィン変換行列を表す
         * @note  このテンプレートから直接特殊化することは想定していない.
         * あくまでmatrixxを実装するための補助テンプレートである.
         * このクラスは不変クラスである.
@@ -32,13 +36,13 @@ namespace math
                > > > > > > > >
        {
        protected:
-//{{{ last_index
+//{{{ size
                /**
-                * @brief components_の最後の添字を得る
+                * @brief components_のサイズを返す
                 */
-               size_t last_index() const
+               size_t size() const
                {
-                       return Dim * (Dim + 1) - 1;
+                       return Dim * (Dim + 1);
                }
 //}}}
                union 
@@ -49,18 +53,27 @@ namespace math
        public:
                // コンストラクタ,デストラクタはデフォルトのものを使う
                // 代入演算子,コピーコンストラクタはデフォルトのものを使う
+//{{{ swap
+               void swap(Derived& rhs) throw()
+               {
+                       using std::swap;
+                       for (size_t i = 0; i < size(); ++i) {
+                               swap(components_[i], rhs.components_[i]);
+                       }
+               }
+//}}}
 //{{{ operator +=
-               Derived& operator+=(const Derived& rhs) const
+               Derived& operator+=(const Derived& rhs)
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] += rhs.components_[i];
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
                }
                
-               Derived& operator+=(float rhs) const
+               Derived& operator+=(float rhs) 
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] += rhs;
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
@@ -68,25 +81,25 @@ namespace math
                
                friend Derived operator+(float rhs1, Derived& rhs2)
                {
-                       float tmp[last_index() + 1];
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       float tmp[size()];
+                       for (size_t i = 0; i < size(); ++i) {
                                tmp[i] = rhs1 + rhs2.components_[i];
                        }
                        return Derived(tmp);
                }
 //}}}
 //{{{ operator -=
-               Derived& operator-=(const Derived& rhs) const
+               Derived& operator-=(const Derived& rhs)
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] -= rhs.components_[i];
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
                }
                
-               Derived& operator-=(float rhs) const
+               Derived& operator-=(float rhs)
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] -= rhs;
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
@@ -94,89 +107,55 @@ namespace math
                
                friend Derived operator-(float rhs1, Derived& rhs2) 
                {
-                       float tmp[last_index() + 1];
-                       for (size_t i = 0; i <= last_index(); ++i) {
-                               retval.components_[i] = rhs1 - rhs2.components_[i];
+                       float tmp[size()];
+                       for (size_t i = 0; i < size(); ++i) {
+                               tmp[i] = rhs1 - rhs2.components_[i];
                        }
                        return Derived(tmp);
                }
 //}}}
 //{{{ operator *=      
-<<<<<<< HEAD
-               /**
-                * @brief 行列の積を計算し,最後の要素が1になるように定数倍する.
-                */
                Derived& operator*=(const Derived& rhs)
                {
-                       Derived retval;
-                       const int SIZE = Dim + 1;
-
-                       // calculate the last element previously
-                       int b = last_index() - Dim;
-                       int c = Dim;
-                       float last_sum = 0;
-                       for (int i = 0; i < SIZE; ++i) {
-                               last_sum += elements_[b + i] * rhs.elements_[c + i * SIZE];
-                       }
-                       retval.elements_[last_index()] = 1;
-
-                       for (int a = last_index() - 1; a >= 0; --a) {
-                               int b = a / SIZE * SIZE;
-                               int c = a % SIZE;
-                               float sum = 0;
-                               for (int i = 0; i < SIZE; ++i) {
-                                       sum += elements_[b + i] * rhs.elements_[c + i * SIZE];
-=======
-               Derived& operator*=(const Derived& rhs) const
-               {
                        Derived M;
-                       for (int i = 0; i < Dim; ++i) {
-                               for (int j = 0; j <= Dim; ++j) {
-                                       for (int k = 0; k <= Dim; ++k) {
-                                               M.elements_[i][j] += (*this)[i][k] * rhs[k][j];
+                       for (size_t i = 0; i < Dim; ++i) {
+                               for (size_t j = 0; j <= Dim; ++j) {
+                                       float sum = 0;
+                                       for (size_t k = 0; k <= Dim; ++k) {
+                                               sum += at(i, k) * rhs.at(k, j);
                                        }
->>>>>>> b5b23396b7e7f34fc19daa767a3cbdad95673d24
+                                       M.elements_[i][j] = sum;
                                }
                        }
                        *this = M;
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
                }
 
-               Derived& operator*=(float rhs) const
+               Derived& operator*=(float rhs) 
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] *= rhs;
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
                }
        
-               Coordinate operator*(const Coordinate& rhs) const
+               Coordinate operator*(const Coordinate& rhs)
                {
                        float arr[Dim];
-                       for (int i = 0; i < Dim; ++i) {
-                               for (int k = 0; k <= Dim; ++k) {
-                                       arr[i] += (*this)[i][k] * rhs[k];
+                       for (size_t i = 0; i < Dim; ++i) {
+                               float sum = 0;
+                               for (size_t k = 0; k <= Dim; ++k) {
+                                       sum += elements_[i][k] * rhs[k];
                                }
+                               arr[i] = sum;
                        }
                        return Coordinate(arr);
                }
-<<<<<<< HEAD
-=======
-
-               friend Derived operator*(float rhs1, Derived& rhs2) 
-               {
-                       float arr[last_index() + 1];
-                       for (size_t i = 0; i <= last_index(); ++i) {
-                               arr[i] *= rhs;
-                       }
-                       return Derived(arr);
-               }
->>>>>>> b5b23396b7e7f34fc19daa767a3cbdad95673d24
 //}}}
 //{{{ operator /=
-               Derived& operator/=(float rhs) const
+               Derived& operator/=(float rhs)
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                components_[i] /= rhs;
                        }
                        return *reinterpret_cast<Derived*>(this);//thisがDerived型であることは保証されている.
@@ -185,7 +164,7 @@ namespace math
 //{{{ operator ==
                bool operator==(const Derived& rhs) const
                {
-                       for (size_t i = 0; i <= last_index(); ++i) {
+                       for (size_t i = 0; i < size(); ++i) {
                                if (std::abs(components_[i] - rhs.components_[i]) > MOF_ERROR_THRESHOLD) return false;
                        }
                        return true;
@@ -196,12 +175,9 @@ namespace math
                 * @note M[i][j]のように参照可能
                 * @note この方法による複数の要素への参照は非効率
                 */
-               const float* const operator [](size_t i) const
+               row_of_matrix<Dim> const operator [](size_t i) const
                {
-                       static float arr[] = {0, 0, 0, 1};
-                       if (i < Dim) return &components_[i * Dim];
-                       return arr;
-                       // TODO index bounds check
+                       return row_of_matrix<Dim>(i, elements_);
                }
 //}}}
 //{{{ operator <<
@@ -236,9 +212,51 @@ namespace math
                 */
                const float at(size_t i, size_t j) const
                {
-                       return elements_[i][j];
+                       if (i < Dim ) return elements_[i][j];
+                       if (j == Dim) return 1;
+                       else return 0;
                }
 //}}}
        };
+//{{{ swap
+       template <size_t Dim, typename Derived, typename Coordinate>
+       void swap
+       (
+               basic_matrix<Dim, Derived, Coordinate>& a,
+               basic_matrix<Dim, Derived, Coordinate>& b
+       ) throw()
+       {
+               a.swap(b);
+       }
+//}}}
+//{{{ row_of_matrix
+       /**
+        * @brief M[i][j]のように行列の要素を取得するための補助クラス
+        */
+       template <size_t Dim>
+       class row_of_matrix : boost::noncopyable
+       {
+               size_t index_;
+               const float (&elements_)[Dim][Dim+1];
+
+       public:
+               row_of_matrix(size_t i, const float (&elements)[Dim][Dim+1])
+                       : index_(i), elements_(elements)
+               {
+               }
+
+               row_of_matrix(const row_of_matrix<Dim>&& rhs)
+                       : index_(rhs.index_), elements_(rhs.elements_)
+               {
+               }
+
+               float operator[](size_t j) const
+               {
+                       if (Dim != index_) return elements_[index_][j];
+                       else if (Dim == j) return 1;
+                       else return 0;
+               }
+       };
+//}}}
 }// namespace math
 }// namespace mof