[an error occurred while processing this directive] [an error occurred while processing this directive][an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] (none) [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive][an error occurred while processing this directive] [an error occurred while processing this directive][an error occurred while processing this directive] [an error occurred while processing this directive][an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive] (none) [an error occurred while processing this directive] [an error occurred while processing this directive] [an error occurred while processing this directive][an error occurred while processing this directive]
 
[an error occurred while processing this directive] [an error occurred while processing this directive]
Skåne Sjælland Linux User Group - http://www.sslug.dk Home   Subscribe   Mail Archive   Forum   Calendar   Search
MhonArc Date: [Date Prev] [Date Index] [Date Next]   Thread: [Date Prev] [Thread Index] [Date Next]   MhonArc
 

Re: [CPROG] lille c++ matrix template lib til openGL brug ??



On Fri, 23 Aug 2002 16:58:35 +0200, Peter Korsgaard wrote:

> http://cvs.48ers.dk/cgi-bin/viewcvs.cgi/inf9/tngcode/inc/

Det ser meget fint ud, men jeg skal bruge lidt mere, idet jeg skal lave
andre ting end openGL i programmet også. Så jeg har også brug for 2x2,
3x3, 3x4, 4x3, 4x4 og lidt til i de små størrelser.

Jeg tænkte på noget i denne stil, og har begyndt på det - men hvis nogen
havde lavet det før mig, var der jo ingen grund til at gøre det igen:

Der mangler en del i svd, lu faktorisering, rotationsmatricer,
og så videre, sorry for den lidt lange posting:

hvis du er interesseret i en mere generel matrix/vektor klasse a' la
denne her, vil jeg gerne slå pjalterne sammen med dig.

hilsen, Marc

-------------------------------------------------------------------
linalg.hpp:

#ifndef LINALG_H
#define LINALG_H


#include <cstddef>
#include <iterator>
#include <algorithm>
#include <math.h>
#include <iostream>
#include <boost/array.hpp>


// forward declaration
// template<class Scalar, std::size_t NR, std::size_t NC> class Mat;


template<class Scalar, std::size_t NR, std::size_t NC> 
class Mat : public boost::array<Scalar, NR*NC> 
{
   public: 
      // type definitions,  notice, inherited from boots::array are all the 
      // usual type definitions
      typedef Scalar                              value_type;
      typedef Scalar*                             iterator;
      typedef const Scalar*                       const_iterator;
      typedef Scalar&                             reference;
      typedef const Scalar&                       const_reference;
      typedef  std::size_t                        size_type;
      typedef  std::ptrdiff_t                     difference_type;
      // now the new type definitons
      typedef  boost::array<std::size_t,2>        index_type;  
      typedef  boost::array<Scalar, NR*NC>        storage_type;
      typedef  boost::array<Scalar, NR*NC>&       storage_reference;
      typedef  const boost::array<Scalar, NR*NC>& const_storage_reference;
      typedef  Mat<Scalar, NR, NC>                container_type;
      typedef  Mat<Scalar, NR, NC>&               container_reference;
      typedef  const Mat<Scalar, NR, NC>&         const_container_reference;
   protected:
      static size_type index(size_type rn, size_type cn) 
      {
	 return cn + rn*NC; // row ordering
	 //return rn + cn*NR; // column ordering
      }
      static size_type index(index_type indarr) 
      {
	 return indarr[1] + indarr[0]*NC; // row ordering
	 //return indarr[0] + indarr[1]*NR; // column ordering
      }
   public:
      // constructors
      Mat() 
      {
	 assign(Scalar());
      }
      //explicit Mat(const Scalar& s = 0) {std::fill_n(begin(), size(), s);};
      explicit Mat(const Scalar& s) 
      {
	 assign(s);
      }
      explicit Mat(const_storage_reference a) 
      {
	 std::copy(a.begin(), a.end(), begin()); 
      }
      Mat(const_container_reference m) 
      {
	 std::copy(m.begin(), m.end(), begin()); 
      }
      // no destructor ~Mat() {};  needed, default is just perfect

      // assignment with type conversion - should I really do this ???
      template <typename Scalar2>
      container_reference operator= (const Scalar2& s) 
      {
	 assign(s); return *this;
      }
      template <typename Scalar2>
      container_reference operator= (const boost::array<Scalar2,NR*NC>& a) 
      {
	 std::copy(a.begin(),a.end(), begin());
	 return *this;
      }
      template <typename Scalar2>
      container_reference operator= (const Mat<Scalar2,NR,NC>& m) 
      {
	 std::copy(m.begin(),m.end(), begin());
	 return *this;
      }
      // size info
      size_type nrows()
      {
	 return NR; 
      }
      size_type ncols()
      {
	 return NC; 
      }
      index_type dim()
      {
	 index_type ind = { { NR, NC } }; 
	 return ind;
      }
      // vector and matrix indexing, offset 0, column or row major, 
      // no index checking
      reference operator[](size_type i)
      {
	 return storage_type::operator[](i);
      }
      const_reference operator[](size_type i) const 
      {
	 return storage_type::operator[](i);
      }
      reference operator[](index_type indarr)
      {
	 return storage_type::operator[](index(indarr));
      }
      const_reference operator[](index_type indarr) const 
      {
	 return storage_type::operator[](index(indarr));
      }

      // vector and matrix indexing, offset 0, column or row major, 
      // no index checking
      reference operator()(size_type i)
      {
	 return storage_type::operator[](i);
      }
      const_reference operator()(size_type i) const 
      {
	 return storage_type::operator[](i);
      }
      reference operator()(size_type rn, size_type cn)
      {
	 return storage_type::operator[](index(rn, cn));
      }
      const_reference operator()(size_type rn, size_type cn) const 
      {
	 return storage_type::operator[](index(rn, cn));
      }
      reference operator()(index_type indarr)
      {
	 return storage_type::operator[](index(indarr));
      }
      const_reference operator()(index_type indarr) const 
      {
	 return storage_type::operator[](index(indarr));
      }

      // vector and matrix indexing, offset 0, column or row major, 
      //index checking
      reference at(size_type i) 
      { 
	 return storage_type::at(i);
      }
      const_reference at(size_type i) const 
      { 
	 return storage_type::at(i);
      }
      reference at(size_type rn, size_type cn)
      {
	 size_type i = index(rn, cn);
	 return storage_type::at(i);
      }
      const_reference at(size_type rn, size_type cn) const 
      {
	 size_type i = index(rn, cn);
	 return storage_type::at(i);
      }
      reference at(index_type indarr)
      {
	 size_type i = index(indarr);
	 return storage_type::at(i);
      }
      const_reference at(index_type indarr) const 
      {
	 size_type i = index(indarr);
	 return storage_type::at(i);
      }
      
      // unary operators
      container_reference operator-() 
      {
	 for(size_type i = 0; i < size(); i++) 
	    elems[i] = -elems[i];
	 return *this;
      }

      // arithmetic matrix operators with type conversion
      template <typename Scalar2>
      container_reference operator+=(const Mat<Scalar2,NR,NC>& m2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] += m2.elems[i];
	 return *this;
      }
      template <typename Scalar2>
      container_reference operator-=(const Mat<Scalar2,NR,NC>& m2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] -= m2.elems[i];
	 return *this;
      }

      // arithmetic scalar operators with type conversion
      template <typename Scalar2>
      container_reference operator+=(const Scalar2& s2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] += s2;
	 return *this;
      }
      template <typename Scalar2>
      container_reference operator-=(const Scalar2& s2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] -= s2;
	 return *this;
      }
      template <typename Scalar2>
      container_reference operator*=(const Scalar2& s2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] *= s2;
	 return *this;
      }
      template <typename Scalar2>
      container_reference operator/=(const Scalar2& s2) 
      {
	 for(size_type i = 0; i < size(); i++) elems[i] /= s2;
	 return *this;
      }

//        // rowvector and colvector operations
//        Mat<Scalar,NR,1>& col(size_type) 
//        {
//  	 Mat <Scalar,NR,1> cvec; 
//  	 return cvec;
//        }
//        Mat<Scalar,NR,1> col(size_type) const 
//        {
//  	 Mat <Scalar,NR,1> cvec; 
//  	 return cvec;
//        }
//        Mat<Scalar,1,NC>& row(size_type) 
//        {
//  	 Mat <Scalar,1,NC> rvec; 
//  	 return rvec;
//        }
//        Mat<Scalar,1,NC> row(size_type) const 
//        {Mat <Scalar,1,NC> rvec; 
//        return rvec;
//        }

      // io read and write operations
      bool read(std::istream& is)  
      {
	 for(size_type i = 0; i < size(); i++) is >> elems[i];
	 return true;
      }
      bool write(std::ostream& os) const 
      {
	 for(size_type i = 0; i < size(); i++) os << elems[i] << " " ;
	 return true;
      }

};  // end of class Mat 


// non-member arithmetic matrix operators with type conversion
template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator+(const Mat<Scalar1,NR,NC>& m1, const Mat<Scalar2,NR,NC>& m2) 
{
   Mat<Scalar1,NR,NC> result = m1;
   return result += m2;
}

template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator+(const Mat<Scalar1,NR,NC>& m, const Scalar2& s) 
{
   Mat<Scalar1,NR,NC> result = m;
   return result += s;
}

template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator+(const Scalar2& s, const Mat<Scalar1,NR,NC>& m) 
{
   Mat<Scalar1,NR,NC> result = m;
   return result += s;
}

template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator*(const Mat<Scalar1,NR,NC>& m, const Scalar2& s) 
{
   Mat<Scalar1,NR,NC> result = m;
   return result *= s;
}

template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator*(const Scalar2& s, const Mat<Scalar1,NR,NC>& m) 
{
   Mat<Scalar1,NR,NC> result = m;
   return result *= s;
}

template <class Scalar1, class Scalar2, std::size_t NR, std::size_t NC>
Mat<Scalar1,NR,NC> 
operator-(const Mat<Scalar1,NR,NC>& m1, const Mat<Scalar2,NR,NC>& m2) 
{
   Mat<Scalar1,NR,NC> result = m1;
   return result -= m2;
}

// multiplication operator
template<class Scalar, std::size_t NR1, std::size_t NCR, std::size_t NC2>
Mat<Scalar,NR1,NC2> operator*(const Mat<Scalar,NR1,NCR>& m1,
			      const Mat<Scalar,NCR,NC2>& m2)
{
   Mat<Scalar,NR1,NC2> result; 
   //  error-prone if not zero-offset used ...
   //  for (Mat<Scalar,NR1,NC2>::size_type rn = 0; rn < NR1; rn++)
   //    for (Mat<Scalar,NR1,NC2>::size_type cn = 0; cn < NC2; cn++)
   //      {
   //        Scalar s = 0.0;
   //	    for (std::size_t crn = 0; crn < NCR; crn++)
   //  	       s += m1(rn,crn)*m2(crn,cn);
   //  	       result(rn,cn) = s;
   //      }

   // not so nice due to clumsy indexing
   //Mat<Scalar,NR1,NC2>::index_type resultInd;
   //Mat<Scalar,NR1,NCR>::index_type m1Ind;
   //Mat<Scalar,NCR,NC2>::index_type m2Ind;
   //for (resultInd[0] = 0; resultInd[0] < NR1; resultInd[0]++)
   //   for (resultInd[1] = 0; resultInd[1] < NC2; resultInd[1]++)
   //   {
   //	 m1Ind[0] = resultInd[0];
   //	 m2Ind[1] = resultInd[1];
   //	 Scalar s = 0.0;
   //	 //for (std::size_t crn = 0; crn < NCR; crn++)
   //	 for (Mat<Scalar,NR1,NC2>::size_type crn = 0; crn < NCR; crn++)
   //	 {
   //	    m1Ind[1] = m2Ind[0] = crn;
   //	    s += m1[m1Ind]*m2[m2Ind];
   //	 }
   //	 result[resultInd] = s;
   //   }

   // last trial: use iterators and advance/distance templates
   for (typename Mat<Scalar,NR1,NC2>::iterator ri = result.begin(); 
   	ri < result.end(); ri++)
   {
      typename Mat<Scalar,NR1,NCR>::const_iterator rowIter = m1.begin();
      // row order
      std::advance(rowIter, NCR*(std::distance(result.begin(),ri)/NC2));
      // col order  
      //std::advance(rowIter, (std::distance(result.begin(),ri)%NR1)); 
      typename Mat<Scalar,NCR,NC2>::const_iterator colIter = m2.begin();
      // row order
      std::advance(colIter, std::distance(result.begin(),ri)%NC2);   
      // col order
      //std::advance(colIter, NCR*(std::distance(result.begin(),ri)/NR1));
      
      Scalar s = 0;
      for (typename Mat<Scalar,NR1,NC2>::size_type crn = 0; crn < NCR; crn++)
      {
	 s += (*rowIter)*(*colIter);
	 std::advance(rowIter,1);    // row ordering
	 std::advance(colIter,NC2);   // row ordering	 
	 //std::advance(rowIter,NR1);  // col ordering
	 //std::advance(colIter,1);   // col ordering	 
      }
      *ri = s;
   }
   return result;
}

// io operators
template<class Scalar, std::size_t NR, std::size_t NC> 
std::ostream& operator<<(std::ostream& os, const Mat<Scalar,NR,NC>& m)
{
   assert(!os.fail());
   m.write(os);
   assert(os.good());
   return os;
}

template<class Scalar, std::size_t NR, std::size_t NC> 
std::istream& operator>>(std::istream& is, Mat<Scalar,NR,NC>& m)
{
   assert(!is.fail());
   m.read(is);
   assert(is.good());
   return is;
}

template<class Scalar, std::size_t NR1, 
         std::size_t NC1, std::size_t NR2, std::size_t NC2>
bool equalDim(Mat<Scalar,NR1,NC1> m1, Mat<Scalar,NR2,NC2> m2)
{
   return (NR1==NR2 && NC1==NC2);
}


//    Mat<Scalar, NR, NC> operator+(const Mat<Scalar,NR,NC>&,
//  					   const Mat<Scalar,NR,NC>&)
//    {Mat<Scalar, NR, NC> mat; return mat;};






// vector definition - partial template specialization typedef does not work...
//template<class Scalar, std::size_t NR> 
//typedef Mat<Scalar,NR,1> ColVec<Scalar,NR>;

//template<class Scalar, std::size_t NC> 
//typedef  Mat<Scalar,1,NC> RowVec<Scalar,NC>;




#endif //LINALG_H


 
Home   Subscribe   Mail Archive   Index   Calendar   Search

 
 
Questions about the web-pages to <www_admin>. Last modified 2005-08-10, 20:09 CEST [an error occurred while processing this directive]
This page is maintained by [an error occurred while processing this directive]MHonArc [an error occurred while processing this directive] # [an error occurred while processing this directive] *