/*
 * LogicalDataType.hh
 * @brief is a complement of ParamData.hh, contain mathematics function between ParamData element
 * @see paramData.hh
 *  Created on: Dec 21, 2012
 *      Author: f.casimir
 */

#ifndef LOGICALDATATYPE_HH_
#define LOGICALDATATYPE_HH_

#include <vector>
#include <cmath>
#include "ParamData.hh"

/**
 * @brief Operator And (&&) between LogicalData
 */
	inline AMDA::Parameters::LogicalData And(AMDA::Parameters::LogicalData a, AMDA::Parameters::LogicalData b) {
		AMDA::Parameters::LogicalData res = AMDA::Parameters::LogicalData::NaN;

		if ( a != AMDA::Parameters::LogicalData::NaN && b != AMDA::Parameters::LogicalData::NaN) {
			res = ((a==AMDA::Parameters::LogicalData::True) && (b==AMDA::Parameters::LogicalData::True))?AMDA::Parameters::LogicalData::True:AMDA::Parameters::LogicalData::False;
		}

		return res;
	}

	/**
	 * @brief Operator Or (||) between LogicalData
	 */
	inline AMDA::Parameters::LogicalData Or(AMDA::Parameters::LogicalData a, AMDA::Parameters::LogicalData b) {
		AMDA::Parameters::LogicalData res = AMDA::Parameters::LogicalData::NaN;

		if ( a != AMDA::Parameters::LogicalData::NaN && b != AMDA::Parameters::LogicalData::NaN) {
			res = (( a == AMDA::Parameters::LogicalData::True) || (b == AMDA::Parameters::LogicalData::True))?AMDA::Parameters::LogicalData::True:AMDA::Parameters::LogicalData::False;
		}

		return res;
	}



	/**
	 * @brief Operator greater_than (>) between LogicalData
	 */
	template <typename Type1, typename Type2>
	AMDA::Parameters::LogicalData greater_than(Type1 a, Type2 b) {
		AMDA::Parameters::LogicalData res = AMDA::Parameters::LogicalData::True;
		if ( std::isfinite(a) && std::isfinite(b)) {
			res = (a>b)?AMDA::Parameters::LogicalData::True:AMDA::Parameters::LogicalData::False;
		} else {
			res = AMDA::Parameters::LogicalData::NaN;
		}
		return res;
	}


	/**
	 * @brief Operator greater_than (>) between table of LogicalData
	 */
	template<typename Type1, typename Type2>
	std::vector<AMDA::Parameters::LogicalData> greater_than(std::vector<Type1> a,
			std::vector<Type2> b) {
		std::vector<AMDA::Parameters::LogicalData> r;
		if (a.size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Logical comparison between two vector of different size"));
		}

		for (unsigned int i = 0; i < a.size(); ++i) {
			r.push_back(greater_than(a[i], b[i]));
		}
		return r;
	}

	/**
	 * @brief Operator lower_than (<) between table of LogicalData
	 */
	template <typename Type1, typename Type2>
	AMDA::Parameters::LogicalData lower_than(Type1 a, Type2 b) {
		AMDA::Parameters::LogicalData res = AMDA::Parameters::LogicalData::True;
		if ( std::isfinite(a) && std::isfinite(b)) {
			res = (a<b)?AMDA::Parameters::LogicalData::True:AMDA::Parameters::LogicalData::False;
		} else {
			res = AMDA::Parameters::LogicalData::NaN;
		}
		return res;
	}

	/**
	 * @brief Operator lower_than (<) between table of LogicalData
	 */
	template<typename Type1, typename Type2>
	std::vector<AMDA::Parameters::LogicalData> lower_than(std::vector<Type1> a,
			std::vector<Type2> b) {
		std::vector<AMDA::Parameters::LogicalData> r;
		if (a.size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Logical comparison between two vector of different size"));
		}

		for (unsigned int i = 0; i < a.size(); ++i) {
			r.push_back(lower_than(a[i], b[i]));
		}
		return r;
	}


	/**
	 * @brief Operator addition between  two  vector
	 * @details
	 * 		r[0] = a[0] + b[0]
	 * 		r[1] = a[1] + b[1]
	 * 		r[2] = a[2] + b[2]
	 * 	@return r
	 */
	template <typename Type1, typename Type2>
	std::vector<Type1> operator +(std::vector<Type1> a, std::vector<Type2> b) {
		std::vector<Type1> r;
		if(a.size() != b.size())  {
			BOOST_THROW_EXCEPTION(AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Add two vector of different size"));
		}

		for(unsigned int i = 0; i < a.size(); ++i) {
			r.push_back((Type1)(a[i] + b[i]));
		}
		return r;
	}


	/**
	 * @brief Operator subtraction between  two  vector
	 * @details
	 * 		r[0] = a[0] - b[0]
	 * 		r[1] = a[1] - b[1]
	 * 		r[2] = a[2] - b[2]
	 * 	@return r
	 */
	template <typename Type1, typename Type2>
	std::vector<Type1> operator -(std::vector<Type1> a, std::vector<Type2> b) {
			std::vector<Type1> r;
		if(a.size() != b.size())  {
			BOOST_THROW_EXCEPTION(AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Subtract two vector of different size"));
		}

		for(unsigned int i = 0; i < a.size(); ++i) {
			r.push_back(a[i] - b[i]);
		}
		return r;
	}



	/**
	 * @brief Operator dot product between  two  vector
	 * @details
	 * 		r = a[0]*b[0] + a[1]*b[1] + a[2]+b[2]
	 * 	@return the scalar r
	 */
	template <typename Type1, typename Type2>
	Type1 operator *(std::vector<Type1> a, std::vector<Type2> b) {
		Type1 r = 0;
		if(a.size() != b.size())  {
			BOOST_THROW_EXCEPTION(AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Multiply two vector of different size"));
		}

		for(unsigned int i = 0; i < a.size(); ++i) {
			r += a[i]*b[i];
		}
		return r;
	}

	/**
	 * @brief Function dot between two vectors can be used for dot product
	 */
        template <typename Type1, typename Type2>
        Type1 dot(std::vector<Type1> a, std::vector<Type2> b) {
                return a*b;
        }

	/**
	 * @brief Operator cross product between  two  vector
	 * @details
	 * 		r[0] = a[1]*b[2] - a[2]*b[1]
	 * 		r[1] = a[2]*b[0] - a[0]*b[2]
	 * 		r[2] = a[0]*b[1] - a[1]*b[0]
	 * 	@return r
	 */
	template <typename Type1, typename Type2>
	std::vector<Type1> operator ^(std::vector<Type1> a, std::vector<Type2> b) {
		std::vector<Type1> r;
		if(a.size() != b.size() || b.size() != 3)  {
			BOOST_THROW_EXCEPTION(AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Vector product on two vector of different size"));
		}

		r.push_back(a[1]*b[2] - a[2]*b[1]);
		r.push_back(a[2]*b[0] - a[0]*b[2]);
		r.push_back(a[0]*b[1] - a[1]*b[0]);

		return r;
	}

	/**
	 * @brief Function cross between two vectors can be used for cross product
	 */
	template <typename Type1, typename Type2>
	std::vector<Type1> cross(std::vector<Type1> a, std::vector<Type2> b) {
		return a^b;
	}

	/**
	 * @brief div each coordinate one to one
	 * @details
	 * 		r[0] = a[0] / b[0]
	 * 		r[1] = a[1] / b[1]
	 * 		r[2] = a[2] /  b[2]
	 */
	template <typename Type1, typename Type2>
	std::vector<Type1> div(std::vector<Type1> a, std::vector<Type2> b) {
		std::vector<Type1> r;
		if(a.size() != b.size())  {
			BOOST_THROW_EXCEPTION(AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Vector div on two vector of different size"));

		}

		for(unsigned int i = 0; i < a.size(); ++i) {
			r.push_back( a[i]/b[i]);
		}

		return r;
	}

	//
	//Idea but more complex: template <typename Type1, typename Type2>
	//auto operator ^(std::vector<Type1> a, std::vector<Type2> b)  -> std::vector<decltype(Type1()*Type2())> {
	//	std::vector<decltype(Type1()*Type2())> r;


//	Idea but more complex: template <typename Type1, typename Type2>
//	double operator ^(Type1 a, Type2 b) {
//		return pow(a, b);
//	}
//	class Integer {
//		public:
//		Integer(int i) : _i(i) {};
//		int _i;
//	};
//
//	inline double operator ^(double a, Integer b) {
//		return pow(a, b._i);
//	}
//
//
//	class Double {
//		public:
//		Double(double i) : _i(i) {};
//		double _i;
//	};
//	inline double operator ^(double a, Double b) {
//		return pow(a, b._i);
//	}

	// vector with scalar operation

	/**
	 * @brief Operator addition between  a  vector and a scalar
	 * @details all member of vector are divided by the scalar
	 */
	template <typename Type>
	std::vector<Type> operator /(std::vector<Type> a, double f) {
		std::vector<Type> r;
		for(typename std::vector<Type>::iterator it = a.begin(); it != a.end(); ++it) {
			r.push_back((Type)(*it / f));
		}
		return r;
	}


	/**
	 * @brief Operator addition between  a  vector and a scalar
	 * @details all member of vector are added with the scalar
	 */
	template <typename Type>
	std::vector<Type> operator +(std::vector<Type> a, double f) {
		std::vector<Type> r;
		for(typename std::vector<Type>::iterator it = a.begin(); it != a.end(); ++it) {
			r.push_back((Type)(*it + f));
		}
		return r;
	}


	/**
	 * @brief Operator subtraction  between a scalar and a  vector
	 * @details all member of vector are subtracted with the scalar
	 */
	template <typename Type>
	std::vector<Type> operator -( double f, std::vector<Type> a) {
		std::vector<Type> r;
		for(typename std::vector<Type>::iterator it = a.begin(); it != a.end(); ++it) {
			r.push_back(f-*it);
		}
		return r;
	}

	/**
	 * @brief Operator addition between a scalar and a  vector
	 * @details all member of vector are added with the scalar
	 */
	template <typename Type>
	std::vector<Type> operator +( double f, std::vector<Type> a) {
		return a + f;
	}

	/**
	 * @brief Operator product between a vector and a scalar
	 * @details all member of vector are added by the scalar
	 */
	template <typename Type1>
	std::vector<Type1> operator *(std::vector<Type1> a, double f) {
		std::vector<Type1> r;
		for(typename std::vector<Type1>::iterator it = a.begin(); it != a.end(); ++it) {
			r.push_back(*it * f);
		}
		return r;
	}

	/**
	 * @brief Operator product between a scalar and a vector
	 * @details all member of vector are added by the scalar
	 */
	template <typename Type>
	std::vector<Type> operator *( double f, std::vector<Type> a) {
		return a * f;
	}



	/**
	 * @return false if val is NaN or inf or -inf
	 */
	template <typename Type>
	inline bool isFinite(Type val) {
		return std::isfinite(val);
	}


	/**
	 * @return false if elem of val is NaN or inf or -inf
	 */
	template <typename Type>
	inline bool isFinite(std::vector<Type> val) {
		bool ret =false;
		auto it = val.begin();
		while(it != val.end() && !ret) {
			ret &= std::isfinite(*it);
			++it;
		}
		return ret;
	}

	/**
	 * @return false if elem of val is NaN or inf or -inf
	 */
	template <typename Type>
	inline bool isFinite(const AMDA::Parameters::Tab2DData<Type>& val) {
		bool ret =false;
		for (int i = 0; i < val.getDim1Size(); ++i)
			for (int j = 0; j < val.getDim2Size(); ++j)
				ret &= std::isfinite(val[i][j]);
		return ret;
	}

	/**
	 * @return the average of std::list
	 */
	template <typename Type>
	inline Type average(std::list<Type> tab, int /*dim1*/, int /*dim2*/) {
		Type mean;
		if(tab.empty()) {
			mean << NotANumber();
		}
		else {
			mean << ElemNull();
			auto it = tab.begin();
			if (it != tab.end()) {
				unsigned int nbNan = 0;
				Type sum;
				sum << NotANumber();
				for (; it != tab.end(); ++it) {
					if(isFinite(*it) && !isNAN(*it)) {
						if (isNAN(sum))
							sum = *it;
						else
							sum = sum + *it;
					}
					else {
						++nbNan;
					}
				}
				if ((tab.size() - nbNan == 0) || isNAN(sum)) {
					mean << NotANumber();
				}
				else {
					mean = sum / (int)(tab.size() - nbNan);
				}
			}
			else {
				mean << NotANumber();
			}
		}
		return mean;
	}

	template <typename Type>
	inline AMDA::Parameters::Tab2DData<Type> average(std::list<AMDA::Parameters::Tab2DData<Type> > tab, int dim1, int dim2) {
		AMDA::Parameters::Tab2DData<Type> mean(dim1,dim2);
		AMDA::Parameters::Tab2DData<unsigned int > nbNan(dim1,dim2);
		AMDA::Parameters::Tab2DData<Type> sum(dim1,dim2);
		if(tab.empty()) {
			mean << NotANumber();
			return mean;
		}
		for (int index1 = 0; index1 < dim1; ++index1) {
			for (int index2 = 0; index2 < dim2; ++index2) {
				nbNan[index1][index2] = 0;
			}
		}
		sum << NotANumber();
		for (auto it = tab.begin(); it != tab.end(); ++it) {
			for(int index1 = 0; index1 < dim1; ++index1) {
				for(int index2 = 0; index2 < dim2; ++index2) {
					Type val = (*it)[index1][index2];
					if(std::isfinite(val) && !isNAN(val)) {
						if (isNAN(sum[index1][index2]))
							sum[index1][index2] = val;
						else
							sum[index1][index2] +=val;
					}
					else {
						++(nbNan[index1][index2]);
					}
				}
			}
		}
		for(int index1 = 0; index1 < dim1; ++index1) {
			for(int index2 = 0; index2 < dim2; ++index2) {
				if (tab.size() - nbNan[index1][index2] == 0) {
					mean[index1][index2] << NotANumber();
				}
				else {
					mean[index1][index2] = sum[index1][index2] / (int)(tab.size() - nbNan[index1][index2]);
				}
			}
		}

		return mean;
	}

	/**
	 * @return the average of std::list
	 */
	template <typename Type>
	inline std::vector<Type> average(std::list<std::vector<Type> > tab, int dim1, int /*dim2*/) {
		std::vector<Type> mean;
		mean.resize(dim1);
		std::vector<unsigned int >nbNan;
		nbNan.resize(dim1);
		std::vector<Type> sum;
		sum.resize(dim1);
		if(tab.empty()) {
			mean << NotANumber();
			return mean;
		}
		for(unsigned int index = 0; index < tab.begin()->size(); ++index) {
			nbNan[index] = 0;
		}
		sum << NotANumber();
		for (auto it = tab.begin(); it != tab.end(); ++it) {
			for(unsigned int index = 0; index < it->size(); ++index) {
				Type val = it->at(index);
				if(std::isfinite(val) && !isNAN(val)) {
					if (isNAN(sum[index]))
						sum[index] = val;
					else
						sum[index] +=val;
				}
				else {
					++(nbNan[index]);
				}
			}
		}
		for(unsigned int index = 0; index < tab.begin()->size(); ++index) {
			if (tab.size() - nbNan[index] == 0) {
				mean[index] << NotANumber();
			}
			else {
				mean[index] = sum[index] / (int)(tab.size() - nbNan[index]);
			}
		}


		return mean;
	}

/**
 *
 */
	template <typename Type>
	inline Type interpolation(const Type& leftVal, const Type& rightVal, double coef, int /*dim1*/, int /*dim2*/) {
		Type val;
		val << NotANumber();
		if (isNAN(leftVal) || isNAN(rightVal) || coef < 0. || coef > 1.)
			return val;
		val = leftVal + (rightVal - leftVal) * coef;
		return val;
        }

/**
 * 
 */

	template <typename Type>
	inline AMDA::Parameters::Tab2DData<Type> interpolation(const AMDA::Parameters::Tab2DData<Type>& leftVal, const AMDA::Parameters::Tab2DData<Type>& rightVal, double coef, int dim1, int dim2) {
		AMDA::Parameters::Tab2DData<Type> val(dim1,dim2);
		val << NotANumber();
		for(int index1 = 0; index1 < dim1; ++index1) {
			for(int index2 = 0; index2 < dim2; ++index2) {
				Type leftValComp = leftVal[index1][index2];
				Type rightValComp = rightVal[index1][index2];
				val[index1][index2] = interpolation(leftValComp, rightValComp, coef, 1, 1);
			}
		}
		return val;
	}

	template <typename Type>
	inline std::vector<Type> interpolation(const std::vector<Type>& leftVal, const std::vector<Type>& rightVal, double coef, int dim1, int /*dim2*/) {
		std::vector<Type> val;
		val.resize(dim1);
		val << NotANumber();
		for(int index = 0; index < dim1; ++index) {
			Type leftValComp = leftVal[index];
			Type rightValComp = rightVal[index];
			val[index] = interpolation(leftValComp, rightValComp, coef, 1, 1);
		}
		return val;
	}

/**
 * @brief Operator + between Tab2D and scalar
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator +(const AMDA::Parameters::Tab2DData<Type1>& a, Type2 b)
{
	AMDA::Parameters::Tab2DData<Type1> r(a);

	if (isNAN(b))
	{
		r << NotANumber();
		return r;
	}

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (!isNAN(r[i][j]))
				r[i][j] += b;
		}
	return r;
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator +(Type1 a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return (b+a);
}

/**
 * @brief Operator - between Tab2D and scalar
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator -(const AMDA::Parameters::Tab2DData<Type1>& a, Type2 b)
{
	return (a+(-b));
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator -(Type1 a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	AMDA::Parameters::Tab2DData<Type2> r(b);

	if (isNAN(a))
	{
		r << NotANumber();
		return r;
	}

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (!isNAN(r[i][j]))
				r[i][j] = a - r[i][j];
		}
	return r;
}

/**
 * @brief Operator * between Tab2D and scalar
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator *(const AMDA::Parameters::Tab2DData<Type1>& a, Type2 b)
{
	AMDA::Parameters::Tab2DData<Type1> r(a);

	if (isNAN(b))
	{
		r << NotANumber();
		return r;
	}

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (!isNAN(r[i][j]))
				r[i][j] *= b;
		}
	return r;
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator *(Type1 a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return (b*a);
}

/**
 * @brief Operator / between Tab2D and scalar
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator /(const AMDA::Parameters::Tab2DData<Type1>& a, Type2 b)
{
	if (b==0)
	{
		AMDA::Parameters::Tab2DData<Type1> r(a);
		r << NotANumber();
		return r;
	}
	return (a*(1./b));
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator /(Type1 a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	AMDA::Parameters::Tab2DData<Type2> r(b);

	if (isNAN(a))
	{
		r << NotANumber();
		return r;
	}

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (isNAN(r[i][j]) || (r[i][j] == 0))
				r[i][j] << NotANumber();
			else
				r[i][j] = a / r[i][j];
		}
	return r;
}

/*
 * @brief Function add between a Tab2D and a vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> add(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b, int mode)
{
	AMDA::Parameters::Tab2DData<Type1> r(a);

	switch (mode)
	{
	case 1 :
		if (r.getDim1Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[i]))
					r[i][j] << NotANumber();
				else
					r[i][j] += b[i];
			}
		break;
	case 2 :
		if (r.getDim2Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[j]))
					r[i][j] << NotANumber();
				else
					r[i][j] += b[j];
			}
		break;
	default:
		BOOST_THROW_EXCEPTION(
			AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}
	return r;
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> add(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b, int mode)
{
	return add(b,a,mode);
}

/**
 * @brief Operator + between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator +(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b)
{
	return add(a,b,1);
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator +(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return add(a,b,1);
}

/**
 * @brief Function sub between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> sub(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b, int mode)
{
	return add(a,(-1)*b,mode);
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> sub(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b, int mode)
{
	AMDA::Parameters::Tab2DData<Type2> r(b);

	switch(mode)
	{
	case 1 :
		if (r.getDim1Size() != a.size()) {
			BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(a[i]))
					r[i][j] << NotANumber();
				else
					r[i][j] = a[i] - r[i][j];
			}
		break;
	case 2 :
		if (r.getDim2Size() != a.size()) {
			BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(a[j]))
					r[i][j] << NotANumber();
				else
					r[i][j] = a[j] - r[i][j];
			}
		break;
	default:
		BOOST_THROW_EXCEPTION(
			AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}
	return r;
}

/**
 * @brief Operator - between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator -(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b)
{
	return sub(a,b,1);
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator -(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return sub(a,b,1);
}

/**
 * @brief Function mult between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> mult(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b, int mode)
{
	AMDA::Parameters::Tab2DData<Type1> r(a);

	switch (mode)
	{
	case 1 :
		if (r.getDim1Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[i]))
					r[i][j] << NotANumber();
				else
					r[i][j] *= b[i];
			}
		break;
	case 2:
		if (r.getDim2Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[j]))
					r[i][j] << NotANumber();
				else
					r[i][j] *= b[j];
			}
		break;
	default:
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}
	return r;
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> mult(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b, int mode)
{
	return mult(b,a,mode);
}

/**
 * @brief Operator * between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator *(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b)
{
	return mult(a,b,1);
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator *(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return mult(b,a,1);
}

/**
 * @brief Function div between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> div(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b, int mode)
{
	AMDA::Parameters::Tab2DData<Type1> r(a);

	switch (mode)
	{
	case 1 :
		if (r.getDim1Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[i]) || (b[i] == 0))
					r[i][j] << NotANumber();
				else
					r[i][j] = r[i][j] / b[i];
			}
		break;
	case 2 :
		if (r.getDim2Size() != b.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(b[j]) || (b[j] == 0))
					r[i][j] << NotANumber();
				else
					r[i][j] = r[i][j] / b[j];
			}
		break;
	default:
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}
	return r;
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> div(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b, int mode)
{
	AMDA::Parameters::Tab2DData<Type2> r(b);

	switch (mode)
	{
	case 1 :
		if (r.getDim1Size() != a.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(a[i]) || (r[i][j] == 0))
					r[i][j] << NotANumber();
				else
					r[i][j] = a[i] / r[i][j];
			}
		break;
	case 2 :
		if (r.getDim2Size() != a.size()) {
			BOOST_THROW_EXCEPTION(
					AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
		}

		for(int i = 0; i < r.getDim1Size(); ++i)
			for(int j = 0; j < r.getDim2Size(); ++j)
			{
				if (isNAN(r[i][j]) || isNAN(a[j]) || (r[i][j] == 0))
					r[i][j] << NotANumber();
				else
					r[i][j] = a[j] / r[i][j];
			}
		break;
	default:
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}
	return r;
}

/**
 * @brief Operator / between Tab2D and vector
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator /(const AMDA::Parameters::Tab2DData<Type1>& a, std::vector<Type2> b)
{
	return div(a,b,1);
}

template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type2> operator /(std::vector<Type1> a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	return div(b,a,1);
}

/**
 * @brief Operator + between two Tab2D
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator +(const AMDA::Parameters::Tab2DData<Type1>& a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	if ((a.getDim1Size() != b.getDim1Size()) || (a.getDim2Size() != b.getDim2Size())) {
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
	}

	AMDA::Parameters::Tab2DData<Type1> r(a);

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (isNAN(a[i][j]) || isNAN(b[i][j]))
				r[i][j] << NotANumber();
			else
				r[i][j] = a[i][j] + b[i][j];
		}

	return r;
}

/**
 * @brief Operator - between two Tab2D
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator -(const AMDA::Parameters::Tab2DData<Type1>& a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	if ((a.getDim1Size() != b.getDim1Size()) || (a.getDim2Size() != b.getDim2Size())) {
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
	}

	AMDA::Parameters::Tab2DData<Type1> r(a);

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (isNAN(a[i][j]) || isNAN(b[i][j]))
				r[i][j] << NotANumber();
			else
				r[i][j] = a[i][j] - b[i][j];
		}
	return r;
}

/**
 * @brief Operator * between two Tab2D
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator *(const AMDA::Parameters::Tab2DData<Type1>& a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	if ((a.getDim1Size() != b.getDim1Size()) || (a.getDim2Size() != b.getDim2Size())) {
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
	}

	AMDA::Parameters::Tab2DData<Type1> r(a);

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (isNAN(a[i][j]) || isNAN(b[i][j]))
				r[i][j] << NotANumber();
			else
				r[i][j] = a[i][j] * b[i][j];
		}
	return r;
}

/**
 * @brief Operator + between two Tab2D
 */
template <typename Type1, typename Type2>
AMDA::Parameters::Tab2DData<Type1> operator /(const AMDA::Parameters::Tab2DData<Type1>& a, const AMDA::Parameters::Tab2DData<Type2>& b)
{
	if ((a.getDim1Size() != b.getDim1Size()) || (a.getDim2Size() != b.getDim2Size())) {
		BOOST_THROW_EXCEPTION(
				AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due to incompatibility in data dimensions"));
	}

	AMDA::Parameters::Tab2DData<Type1> r(a);

	for(int i = 0; i < r.getDim1Size(); ++i)
		for(int j = 0; j < r.getDim2Size(); ++j)
		{
			if (isNAN(a[i][j]) || isNAN(b[i][j]) || (b[i][j] == 0))
				r[i][j] << NotANumber();
			else
				r[i][j] = a[i][j] / b[i][j];
		}
	return r;
}

/**
 * @brief Function total on all components of a Tab2D => result is a scalar
 */
template <typename Type>
Type total(const AMDA::Parameters::Tab2DData<Type>& a)
{
	Type sum;
	sum << NotANumber();

	for(int i = 0; i < a.getDim1Size(); ++i)
		for(int j = 0; j < a.getDim2Size(); ++j)
		{
			if (isFinite(a[i][j]))
			{
				if (isNAN(sum))
					sum = a[i][j];
				else
					sum += a[i][j];
			}
		}
	return sum;
}

/**
 * @brief Function total on lines or columns of a Tab2D => result is a vector
 * mode = 1 => sum on lines
 * mode = 2 => sum on columns
 */
template <typename Type>
std::vector<Type> total(const AMDA::Parameters::Tab2DData<Type>& a,int mode)
{
	std::vector<Type> sum;

	switch (mode)
	{
	case 1 :
		for(int j = 0; j < a.getDim2Size(); ++j)
		{
			Type elem;
			elem << NotANumber();
			for(int i = 0; i < a.getDim1Size(); ++i)
			{
				if (isFinite(a[i][j]))
				{
					if (isNAN(elem))
						elem = a[i][j];
					else
						elem += a[i][j];
				}
			}
			sum.push_back(elem);
		}
		break;
	case 2 :
		for(int i = 0; i < a.getDim1Size(); ++i)
		{
			Type elem;
			elem << NotANumber();
			for(int j = 0; j < a.getDim2Size(); ++j)
			{
				if (isFinite(a[i][j]))
				{
					if (isNAN(elem))
						elem = a[i][j];
					else
						elem += a[i][j];
				}
			}
			sum.push_back(elem);
		}
		break;
	default:
		BOOST_THROW_EXCEPTION(
			AMDA::AMDA_exception() << AMDA::errno_code(AMDA_OPER_NOT_ALLOWED) << AMDA::ex_msg("Operation not allowed due a not implemented mode"));
	}

	return sum;
}



#endif /* LOGICALDATATYPE_HH_ */