//Matrix class
//Matrix class for storing training and test data
//feature extraction
#pragma once
#ifndef MATRIX_H
#define MATRIX_H

#include <vector>
#include <numeric>
#include <limits>
#include <iterator>
#include <algorithm>

using std::vector;
using std::size_t;


//Matrix class
//Simple class which represents 2 dimensional matrix.
//Stores the data continuously in a vector (row major ordering)
//Examle of usage:
//      //Will create matrix with 2 rows and 4 columns
//      Matrix<double> mat({0, 1, 2, 3,
//                      4, 5, 6 ,7}, 2, 4)
//
//Element access throgh operator():
//      double x = mat(1,2) 
//
//Acess to the row also through operator() 
//      vector<double> row = mat(1) 
//Will store in row {4,5,6,7}
//TODO:
//There is a problem with fast access to the row
template<class T>
class Matrix {
private:
    //matrix shape, number of rows and columns
    size_t  _n_rows, _n_cols;
	//vector to store the data, (row-major ordering)
	vector<T> _data;
public:
	const vector<T>& data() const {
		return _data;
	}

    //define default c'tor?
    Matrix() : _n_rows(0), _n_cols(0) {}

    //init matrix with defined shape
    Matrix(size_t  n_rows, size_t  n_cols);
    //init matrix with defined shape and value
    Matrix(size_t  n_rows, size_t  n_cols, T init_val);
    //create matrix from std::vector with defined shape, (row-major ordering)
    Matrix(vector<T> data, size_t  n_rows, size_t  n_cols);
    Matrix(vector<T> data, size_t  n_cols);

    //create matrix from raw data (row-major ordering)
    //with defined shape
    Matrix(T *data, size_t  n_rows, size_t  n_cols);

    //n_rows getter
    const std::size_t& n_rows() const { return _n_rows; }

    //n_rows getter
    const std::size_t& n_cols() const { return _n_cols; }

    //element setter
    T& operator()(size_t  i, size_t  j);

    //element getter
    const T& operator()(size_t  i, size_t  j) const;

    //row getter
    const typename vector<T> operator()(size_t  i) const;

    //row setter
    Matrix& set_row(const vector<T>& row, size_t  row_index);
    
    Matrix<T>& add_row(const vector<T>& row);

    //return the maximum along a given axis
    vector<T> max(size_t  axis = 0) const;

    //return the maximum along a given axis
    vector<T> min(size_t  axis = 0) const;
};

template <class T>
inline
Matrix<T>::Matrix(size_t  n_rows, size_t  n_cols) :
_n_cols(n_cols), _n_rows(n_rows), _data(_n_rows*_n_cols) {
    ;
}


template <class T>
inline
Matrix<T>::Matrix(size_t  n_rows, size_t  n_cols, T init_val) :
_n_cols(n_cols), _n_rows(n_rows), _data(_n_rows*_n_cols, init_val) {
    ;
}

template <class T>
inline
Matrix<T>::Matrix(vector<T> data, size_t  n_rows, size_t  n_cols) :
_data(data), _n_cols(n_cols), _n_rows(n_rows) {
    if(_n_cols*_n_rows != data.size())
        throw std::length_error("wrong data shape");
}

template <class T>
inline
Matrix<T>::Matrix(vector<T> data, size_t n_cols) :
_data(data), _n_cols(n_cols), _n_rows(data.size() / n_cols) {
    if(_n_cols*_n_rows != data.size())
        throw std::length_error("wrong data shape");
}


template <class T>
inline
Matrix<T>::Matrix(T *data, size_t  n_rows, size_t  n_cols) :
_n_cols(n_cols), _n_rows(n_rows) {
    this->_data.assign(data, data + _n_rows*_n_cols);
}


template <class T>
inline T&
Matrix<T>::operator()(size_t  i, size_t  j) {
    return this->_data[i*this->_n_cols + j];
}

template <class T>
inline const T&
Matrix<T>::operator()(size_t  i, size_t  j) const {
    return this->_data[i*this->_n_cols + j];
}


//Need to do something with this, probably someting like this:
//
//template <class T>
//struct Row{
//    T *data;
//    size_t  len;
//};
//
//template <class T>
//inline const Row
//Matrix<T>::operator()(size_t  i) const
//{
//    Row row;
//    row.data = &this->_data[i*this->_n_cols];
//    row.len = this->_n_cols();
//    return row;
//}
template <class T>
inline const
typename vector<T> //should be const?
Matrix<T>::operator()(size_t  i) const {
    typename vector<T>::const_iterator begin = this->_data.begin() + i*this->_n_cols;
    vector<T> vect(begin, begin + this->_n_cols);
    return vect;
}

template <class T>
Matrix<T>&
Matrix<T>::set_row(const vector<T>& row, size_t  row_index) {
    if(row.size() != _n_cols)
        throw std::out_of_range("number of elements is wrong");
    for(size_t i = 0; i < n_cols(); i++) {
        (*this)(row_index, i) = row[i];
    }
    return (*this);
}

template <class T>
Matrix<T>&
Matrix<T>::add_row(const vector<T>& row) {
    if(row.size() != _n_cols)
        throw std::out_of_range("number of elements is wrong");
    _n_rows++;
    _data.resize(_data.size() + _n_cols);
    for(size_t i = 0; i < row.size(); i++)
        (*this)(_n_rows - 1, i) = row[i];
    return *this;
}

template <class T>
vector<T>
Matrix<T>::max(size_t  axis) const {
    if(axis == 1) {
        vector<T> result(this->_n_rows);
        for(size_t i = 0; i < this->_n_rows; i++) {
            //find maximum at current row
            result[i] = *std::max_element(
                this->_data.begin() + i*this->_n_cols,
                this->_data.begin() + (i + 1)*this->_n_cols);
        }
        return result;
    }
    else {
        vector<T> result(this->_n_cols);
        for(size_t i = 0; i < this->_n_cols; i++) {
            //iterate through column and find max (at current column)
            T tmax = std::numeric_limits<T>::min();
            for(size_t j = 0; j < this->_n_rows; j++) {
                if(tmax < (*this)(j, i))
                    tmax = (*this)(j, i);
            }
            result[i] = tmax;
        }
        return result;
    }

}


template <class T>
vector<T>
Matrix<T>::min(size_t  axis) const {
    if(axis == 1) {
        vector<T> result(this->_n_rows);
        for(size_t i = 0; i < this->_n_rows; i++) {
            //find maximum at current row
            result[i] = *std::min_element(
                this->_data.begin() + i*this->_n_cols,
                this->_data.begin() + (i + 1)*this->_n_cols);
        }
        return result;
    }
    else {
        vector<T> result(this->_n_cols);
        for(size_t i = 0; i < this->_n_cols; i++) {
            //iterate through column and find min (at current column)
            T tmin = std::numeric_limits<T>::max();
            for(size_t j = 0; j < this->_n_rows; j++) {
                if(tmin > (*this)(j, i))
                    tmin = (*this)(j, i);
            }
            result[i] = tmin;
        }
        return result;
    }

}


#endif

