/*********************************************************************
 * This file is part of the cpplibs suite.
 *
 * Copyright (C) 2001 Topi Mäenpää and Jaakko Viertola
 * All rights reserved.
 *
 * This program is free software. You can redistribute and/or modify
 * it under the terms of the free software licence found in the
 * accompanying file "COPYING". The licence terms must always be
 * redistributed with this source file. The above copyright notice
 * must be reproduced in all modified and unmodified copies of this
 * source file.
 *
 * $Revision: 1.12 $
 *********************************************************************/

#ifndef _MATRIXUTILS_H
#define _MATRIXUTILS_H

#include <Matrix.h>

namespace util
{
	/**
	 * Methods for non-mathematics related matrix operations.
	 **/
	class MatrixUtils
	{
	public:
		/**
		 * Possible types of matrix extensions.
		 *
		 * <ul>
		 * <li>EXTEND_ZEROPAD - pad with zeros
		 * <li>EXTEND_REPLICATE - replicate the value on the border
		 * <li>EXTEND_SYMMETRIC - symmetrically mirror boundary values
		 * <li>EXTEND_PERIODIC - treat the matrix as periodic
		 * </ul>
		 **/
		enum ExtendType { EXTEND_ZEROPAD,
											EXTEND_REPLICATE,
											EXTEND_SYMMETRIC,
											EXTEND_PERIODIC };


		/**
		 * A constant for the top of a matrix (1).
		 **/
		static int TOP;
		/**
		 * A constant for the bottom of a matrix (2).
		 **/
		static int BOTTOM;
		/**
		 * A constant for the left side of a matrix (4).
		 **/
		static int LEFT;
		/**
		 * A constant for the right side of a matrix (8).
		 **/
		static int RIGHT;

		/**
		 * Extend a matrix by adding data to its sides.
		 *
		 * @param mat the matrix to be extended
		 * @param amount the number of items to add to each side
		 * @param type the extension type
		 * @param sides a bit mask that tells which sides of the matrix
		 * are extended. Use the constant values TOP, BOTTOM, LEFT and
		 * RIGHT. TOP+BOTTOM = top and bottom and so on.
		 **/
		template <class T> static Matrix<T> extend(const Matrix<T>& mat,
																							 int amount,
																							 ExtendType type = EXTEND_ZEROPAD,
																							 int sides = -1);

		/**
		 * Turn a matrix around so that its top left value is moved to the
		 * bottom right corner. The effect is the same as calling flipud()
		 * and fliplr() sequentially. The matrix contents are modified in
		 * place.
		 **/
		template <class T> static void turnAround(Matrix<T>& mat);

		/**
		 * Split a matrix into sub-matrices. The maximum number of
		 * sub-matrices is taken starting from the given upper left
		 * corner. The sub-matrices are stored in a list in horizontal
		 * raster-scan order.
		 *
		 * @param mat the matrix to be split
		 * @param rows the height of the sub-matrices
		 * @param cols the width of the sub-matrices
		 * @param startr the row of the start position (upper left corner)
		 * @param startc the column of the start position (upper left corner)
		 **/
		template <class T> static List<Matrix<T> > split(const Matrix<T>& mat,
																										 int rows, int cols, int startr = 0, int startc = 0);

		/**
		 * Compare matrix elements to a value. All matrix items that meet
		 * the comparison requirement are either set to one (setToOne ==
		 * true) or left intact (setToOne == false). Items that do not
		 * meet the requirement are set to zero. The size of the returned
		 * matrix is the same as the size of the input matrix. Examples:
		 *
		 * <pre>
		 * #include &lt;functional&gt;
		 * ...
		 * Matrix&lt;int&gt; test(3,3,
		 *                  1,2,3,
		 *                  4,5,6,
		 *                  7,8,9);
		 * //Retain all items that are larger than 5 (set others to 0)
		 * Matrix&lt;int&gt; large(MatrixUtils::compare&lt;std::greater&lt;int&gt; &gt;(test, 5, false);
		 *
		 * //Set six to one and others to zero
		 * Matrix&lt;int&gt; six(MatrixUtils::compare&lt;std::equal_to&lt;int&gt; &gt;(test, 6));
		 * </pre>
		 **/
		template <class comparator, class T> static Matrix<T> compare(const Matrix<T>& mat,
																																	T value,
																																	bool setToOne = true);

		/**
		 * Collect the coordinates of all items in a matrix meeting a
		 * comparison requirement. The values in <i>mat</i> are compared
		 * to <i>value</i> with the comparator given as a template
		 * parameter. If a match is found, the corresponding row and
		 * column indices are placed into the two lists supplied as
		 * parameters. Example:
		 *
		 * <pre>
		 * #include &lt;functional&gt;
		 * ...
		 * Matrix&lt;int&gt; test(3,3,
		 *                  1,2,3,
		 *                  4,5,6,
		 *                  7,8,9);
		 * List<int> rows, cols;
		 * //Collect the coordinates of all items that are smaller than or equal to three
		 * MatrixUtils::findCoordinates&lt;std::less_equal&lt;int&gt; &gt;(test, 3, rows, cols);
		 * //rows = {0, 0}, cols = {0, 1}
		 * </pre>
		 **/
		template <class comparator, class T> static void findCoordinates(const Matrix<T>& mat,
																																		 T value,
																																		 List<int>& rows,
																																		 List<int>& columns);


		/**
		 * Copy a matrix on top of another. The contents of the source
		 * matrix replace the values in the target matrix, starting at the
		 * given upper left row and column coordinates. If the source
		 * matrix is too large to fit inside the target, the parts that
		 * fall outside the boundaries of the target are ignored.
		 *
		 * @param source the matrix to be copied
		 * @param target the matrix on which source is copied
		 * @param targetR the start row at target
		 * @param targetC the start column at target
		 **/
		template <class T> static void copy(const Matrix<T>& source,
																				Matrix<T>& target,
																				int targetR, int targetC);
	};

	template <class comparator, class T> Matrix<T> MatrixUtils::compare(const Matrix<T>& mat,
																																			T value,
																																			bool setToOne)
	{
		comparator comp;
		int r = mat.getRows(), c = mat.getColumns();
		Matrix<T> result(r, c, false);
		T* tData = result.getData();
		const T* sData = mat.getData();
		for (int i=r*c; i--; tData++, sData++)
			{
				if (!comp(*sData, value))
					*tData = T(0);
				else if (setToOne)
					*tData = T(1);
				else
					*tData = *sData;
			}
		return result;
	}

	template <class comparator, class T> void MatrixUtils::findCoordinates(const Matrix<T>& mat,
																																				 T value,
																																				 List<int>& rows, List<int>& cols)
	{
		comparator comp;
		int r,c;
		AllItems(r, c, mat)
			if (comp(mat(r,c), value))
				{
					rows += r;
					cols += c;
				}
	}

	template <class T> Matrix<T> MatrixUtils::extend(const Matrix<T>& mat,
																									 int amount,
																									 ExtendType type,
																									 int sides)
	{
		int oldR = mat.getRows(), oldC = mat.getColumns();
		int doubleAmount = amount<<1;
		int leftCols = 0, rightCols = 0, topRows = 0, bottomRows = 0;

		if (sides & TOP)
			topRows = amount;
		if (sides & BOTTOM)
			bottomRows = amount;
		if (sides & LEFT)
			leftCols = amount;
		if (sides & RIGHT)
			rightCols = amount;

		int newR = oldR+topRows+bottomRows, newC = oldC+leftCols+rightCols;

		//Make all elements zero if extension type is zero padding
		Matrix<T> result(newR, newC, type == EXTEND_ZEROPAD);

		//Copy inner contents
		const T* ps = mat.getData();
		T *p1, *p2 = result.getData()+topRows*newC+leftCols;
		for (int i=oldR; i--; p2+=leftCols+rightCols)
			for (int j=oldC; j--; ps++,p2++)
				*p2 = *ps;
				
		switch (type)
			{
			case EXTEND_ZEROPAD:
				break;

			case EXTEND_SYMMETRIC:
				if (topRows)
					{
						p2 = result.getData()+leftCols;
						p1 = p2 + (doubleAmount-1) * newC;
						for (int i=amount; i--; p1-=(newC<<1)-(rightCols+leftCols),p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p2 = *p1;
					}

				if (bottomRows)
					{
						p2 = result.getData()+(newR-doubleAmount)*newC+leftCols;
						p1 = p2 + (doubleAmount-1)*newC;
						for (int i=amount; i--; p1-=(newC<<1)-(rightCols+leftCols),p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p1 = *p2;
					}

				p2 = result.getData();
				p1 = p2 + doubleAmount-1;
				for (int i=newR; i--; p1+=(rightCols+leftCols)<<1)
					{
						if (leftCols)
							{
								for (int j=leftCols; j--; p1--,p2++)
									*p2 = *p1;
								p2 += oldC;
								p1 += oldC;
							}
						else
							{
								p2 += oldC;
								p1 = p2-1;
							}
						for (int j=rightCols; j--; p1--,p2++)
							*p2 = *p1;
					}

				break;
				
			case EXTEND_PERIODIC:
				if (topRows)
					{
						p2 = result.getData()+leftCols;
						p1 = p2 + oldR*newC;
						for (int i=amount;i--; p1+=rightCols+leftCols, p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p2 = *p1;
					}

				if (bottomRows)
					{
						p1 = result.getData()+leftCols+(amount)*newC;
						p2 = p1 + oldR*newC;
						for (int i=amount; i--; p1+=rightCols+leftCols, p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p2 = *p1;
					}

				p2 = result.getData();
				for (int i=newR; i--; p1+=leftCols)
					{
						p1 = p2+oldC;

						for (int j=leftCols; j--; p1++, p2++)
							*p2 = *p1;

						p1 = p2;
						p2 += oldC;

						for (int j=rightCols; j--; p1++, p2++)
							*p2 = *p1;
					}

				break;

			case EXTEND_REPLICATE:
				if (topRows)
					{
						p2 = result.getData()+leftCols;
						p1 = p2 + amount*newC;
						for (int i=amount;i--; p1-=oldC, p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p2 = *p1;
					}

				if (bottomRows)
					{
						p1 = result.getData()+(newR-amount-1)*newC+leftCols;
						p2 = p1 + newC;
						for (int i=amount; i--; p1-=oldC, p2+=rightCols+leftCols)
							for (int j=oldC; j--; p1++, p2++)
								*p2 = *p1;
					}

				p2 = result.getData();
				for (int i=newR; i--; p1+=leftCols)
					{
						p1 = p2+leftCols;

						for (int j=leftCols; j--; p2++)
							*p2 = *p1;

						p2 += oldC;
						p1 += oldC-1;

						for (int j=rightCols; j--; p2++)
							*p2 = *p1;
					}
				break;
				
			default:
				break;
			}

		return result;
	}

	template <class T> void MatrixUtils::turnAround(Matrix<T>& mat)
	{
		int rows = mat.getRows(), cols = mat.getColumns(), size = rows*cols;
		T* data1 = mat.getData(), *data2 = data1 + size-1;
		T tmp;
		for (int i=size>>1; i--; data1++, data2--)
			{
				tmp = *data1;
				*data1 = *data2;
				*data2 = tmp;
			}
	}

	template <class T> List<Matrix<T> > MatrixUtils::split(const Matrix<T>& mat,
																												 int rows, int cols, int startr, int startc)
	{
		int r = mat.getRows(), c = mat.getColumns();
		int nr = (r-startr) / rows;
		int nc = (c-startc) / cols;
		
		List<Matrix<T> > result(nr*nc);
		for (int r=0;r<nr;r++)
			for (int c=0;c<nc;c++)
				result += mat(r*rows+startr,c*cols+startc,rows,cols);

		return result;
	}

	template <class T> void MatrixUtils::copy(const Matrix<T>& source,
																						Matrix<T>& target,
																						int targetR, int targetC)
	{
		int sRows = source.getRows(), sCols = source.getColumns();
		int tRows = target.getRows(), tCols = target.getColumns();
		int rowLength = targetC + sCols > tCols ? tCols - targetC : sCols;
		int rowCount = targetR + sRows > tRows ? tRows - targetR : sRows;

		const T* sData = source.getData();
		T* tData = target.getData() + targetR*tCols + targetC;

		if (rowCount > 0 && rowLength > 0)
			for (int i=rowCount; i--; sData += sCols, tData += tCols)
				Util::copyArray(sData, tData, rowLength);
	}
}

#endif