matlab.hpp

#ifndef MATLAB_HPP_
#define MATLAB_HPP_

#include <iostream>
#include <boost/numeric/ublas/blas.hpp>
#include <boost/numeric/ublas/symmetric.hpp>
#include "jmath/indirectArray.hpp"

namespace ublas = boost::numeric::ublas;

namespace jafar {
  namespace jmath {


    //   CLASS MATLAB

    class MATLAB {

      public:
        static bool fullPrec;
        static const size_t prec;
        std::string str;
        friend std::ostream & operator <<(std::ostream & os, const MATLAB & m) {
          return os << m.str;
        }

        template<typename bubTemplateVector>
        void initFromBubVector(const bubTemplateVector& v1) {
          std::ostringstream os;
          os.precision(prec);
          //          os.setf(std::ios::fixed,std::ios::floatfield);
          os.setf(std::ios::fixed,std::ios::adjustfield);
          os << "[ ";
          for (size_t i = 0; i < v1.size(); ++i) {
            {
              os << " ";
              double a = v1(i);
              os << a;
            }
            if (v1.size() != i + 1) {
              os << ", ";
            }
          }
          os << " ]';";
          str = os.str();
        }

        template<typename bubTemplateMatrix>
        void initFromBubMatrix(const bubTemplateMatrix& m1) {
          fullPrec = false;
          std::ostringstream os;
          os << "...\n[ ";
          std::ostringstream ostmp;
          ostmp.precision(prec);
//          ostmp.setf(std::ios::fixed,std::ios::floatfield);
          ostmp.setf(std::ios::fixed,std::ios::adjustfield);
          for (size_t i = 0; i < m1.size1(); ++i) {
            for (size_t j = 0; j < m1.size2(); ++j) {
              if (m1(i, j) < 0)
                ostmp << "-";
              else
                ostmp << " ";
              if (fullPrec || fabs(m1(i, j)) > 1e-6)
                ostmp << fabs(m1(i, j));
              else {
                ostmp << "0";
              }
              if (m1.size2() != j + 1) {
                ostmp << ",";
                const int size = ostmp.str().length();
                for (size_t i = size; i < prec+5; ++i)
                  ostmp << " ";
              }
              os << ostmp.str();
              ostmp.str("");
            }
            if (m1.size1() != i + 1) {
              os << " ;" << std::endl << "  ";
            }
            else {
              os << " ];";
            }
          }
          str = os.str();
        }

        template<typename bubTemplateIndex>
        void initFromBubIndex(const bubTemplateIndex & i1) {
          std::ostringstream os;
          os << "[ ";
          for (size_t i = 0; i < i1.size(); ++i) {
            os << (i1(i) + 1);
            if (i1.size() != i + 1) { // Matlab indexes from 1 to size()
              os << ", ";
            }
          }
          os << " ];";
          str = os.str();
        }


        // Specialized templates for vectors

        template<typename T>
        MATLAB(const ublas::vector<T>& v1) {
          initFromBubVector(v1);
        }

        template<typename T, std::size_t N>
        MATLAB(const ublas::bounded_vector<T, N>& v1) {
          initFromBubVector(v1);
        }


        // Specialized templates for indirect arrays, ranges and slices

        MATLAB(const ublas::indirect_array<>& i1) {
          initFromBubIndex(i1);
        }

        MATLAB(const ublas::range& r1) {
          initFromBubIndex(jafar::jmath::ublasExtra::ia_set(r1));
        }

        MATLAB(const ublas::slice& s1) {
          initFromBubIndex(jafar::jmath::ublasExtra::ia_set(s1));
        }


        // Specialized templates for vector proxies

        template<typename T>
        MATLAB(const ublas::vector_indirect<T>& v1) {
          initFromBubVector(v1);
        }

        template<typename T>
        MATLAB(const ublas::vector_range<T>& v1) {
          initFromBubVector(v1);
        }

        template<typename T>
        MATLAB(const ublas::vector_slice<T>& v1) {
          initFromBubVector(v1);
        }


        // Non-specialized templates will be interpreted as matrices.
        // These include all kinds of 2-dimensional arrays:

        template<class M>
        MATLAB(const M& m1) {
          initFromBubMatrix(m1);
        }

    };

  }
}

#endif /* MATLAB_HPP_ */