/* $Id */

/**
 * @file
 * Implements a Fixed-point data type.
 * This files defines a new data type
 * "FixedT<Tstorage, Tdec_bits>", which can be used to compute fractional
 * numbers without having to worry about network stability (because
 * there is no real rounding taking place.
 * Use "Fixed" for the defaults values (48 bits integral part, 16 bits
 * fractional part) or "FixedT<typename, int>" for flexible sizes. "int"
 * number of bits are taken from a variable of "typename" for the fractional
 * representation, the rest stays available for the integral part.
 * @warning More than 16 bits fractional should not be used when working
 * with large integers.
 * @warning More than 31 fractional bits are not supported, and will
 * trigger a compiler warning with appropriate settings.
 * @warning There are no warnings for buffer overflows. Those are very likely
 * to happen with divisions, as those first shift the numerator by the number
 * of fractional bits. Use multiplication with the reciprocal whereever
 * it is possible.
 * @warning Please use signed variables for storage at the moment ONLY
 *
 * Variables can be initialized:
 * 	- with integers: "FixedT<int, 9> a = -5"
 * 	- with fractionals: "FixedT<int64, 11> b(7, 2);" where the first number is
 * the numerator and the second is the denominator.
 * 	- not at all: "Fixed c;", which will automatically set c to zero.
 *
 * @example fixed_ex.cpp How to initialize Fixed-point data types
 */

/* Included to overload the stream operator << */
#include <iostream>
/* Needed for some ottd-specific data types */
#include "stdafx.h"

/**
 * Base Class for storing fixed-point data types.
 * Fixed-point data types are stored in a single integer-type variable
 * (Tstorage) of which a certain number of bits (Tdec_bits) is used to
 * represent the fractional part of the number.
 */
template <typename Tstorage, int Tdec_bits>
class FixedRawT {
public:
	/** The number of bits that represent the fractional */
	static const int dec_bits = Tdec_bits;

	/** The storage of the number itself */
	Tstorage m_data;

	/**
	 * Basic constructor for integer arguments.
	 * @param value The raw value to which we want to initialize
	 */
	template <typename T>
	FixedRawT(T value) : m_data(value) {}

	/**
	 * Make sure that a Fixed variable is always inizialized to zero by default.
	 * We need this ctor so that we can declare variable without
	 * the need to initialize it explicitly
	 */
	FixedRawT() : m_data (0) {}

	/**
	 * Returns the biggest integral number we can represent
	 */
	int64 FIXED_MAX() const { return m_maximum; }

	/**
	 * Returns the biggest integral number we can represent
	 */
	int64 FIXED_MIN() const { return m_minimum; }

private:
	/** The largest number we can hold */
	static const int64 m_maximum = +(1ULL << ((sizeof(Tstorage) * 8) - Tdec_bits - 1)) - 1;

	/** The smallest number we can hold */
	static const int64 m_minimum = -(1ULL << ((sizeof(Tstorage) * 8) - Tdec_bits - 1)) - 0;
};

/* forward-declare some structs */
template <typename Tstorage, int Tdec_bits> class FixedT;
template <typename T> struct FixedHelperT;

/**
 * Specialization of FixedHelperT. Used to decimal-align two variables of Fixed type.
 */
template <typename Tstorage, int Tdec_bits> struct FixedHelperT<FixedRawT<Tstorage, Tdec_bits> >
{
	/** The number of bits used for the fraction */
	static const int dec_bits = Tdec_bits;

	/**
	 * Gives the raw data of a FixedRaw
	 * @param t The number to return
	 */
	static int64 Raw(const FixedRawT<Tstorage, Tdec_bits>& t) {return (int64)t.m_data;}

	/**
	 * Returns the data from FixedRawT aligned so that it is aligned to a given number
	 * @param t    The number to be aligned
	 * @param bits The number of bits in the fractional part to be aligned to
	 */
	static int64 Raw(const FixedRawT<Tstorage, Tdec_bits>& t, int bits)
	{
		return (((int64)t.m_data) << (bits > Tdec_bits ? bits - Tdec_bits : 0)) >> (bits < Tdec_bits ? Tdec_bits - bits : 0);
	}
};

/**
 * Specialization of FixedHelperT. Used to decimal-align two variables of Fixed type.
 */
template <typename Tstorage, int Tdec_bits> struct FixedHelperT<FixedT<Tstorage, Tdec_bits> > : public FixedHelperT<FixedRawT<Tstorage, Tdec_bits> > {};

/**
 * General implementation of FixedHelperT. This tempate makes sure that
 * some number supplied is properly aligned at the decimal
 */
template <typename T> struct FixedHelperT
{
	/** This version of FixedHelperT is only used for full integers, so assume the number of dec_bits to be zero */
	static const int dec_bits = 0;
	/** Converts a full integer to one that has a given number of fractional bits, assumes
	 * the number of fractional bits to be zero.
	 * @param t the number to be converted
	 */
	static int64 Raw(const T& t) {return (int64)t;}

	/**
	 * Converts a full integer to one that has a given number of fractional bits
	 * @param bits the number of bits for the fraction
	 * @param t    the number to be converted
	 */
	static int64 Raw(const T& t, int bits) {return ((int64)t) << bits;}
};


/**
 * A class that defines a fixed-point data type, which a variable length and precision.
 * The data type that is defined is a fixed-point variable that has and number of Tdec_bits
 * bits to represent the fractional part, the remaining bits of Tstorage are then used for
 * the integer part of the number. This means, we always need to make a trade-off between
 * the precision we want (higher number in Tdec_bits) and the range (smaller number of
 * Tdec_bits). This class here only defines the constructors and the operators, everything
 * else is happning in helper classes.
 * @note for usage of fprintf and similar, explicit casts are needed (double or int64).
 */
template <typename Tstorage, int Tdec_bits> class FixedT : public FixedRawT<Tstorage, Tdec_bits> {
/* private block up here, because we need the typedef later on */
private:
	/** We shortcut the underlying data type to "Raw", to save typing */
	typedef FixedRawT<Tstorage, Tdec_bits> Raw;

public:
	/**
	 * Ctor for assignment with other, non floating, variable types
	 * @param value The Value we should assign to the Fixed
	 */
	template <typename T> FixedT(T value) : Raw(FixedHelperT<T>::Raw(value, Tdec_bits)) {}

	/** As a constructor without initializing, just use the one that Raw uses. Sets a variable to zero */
	FixedT() : Raw() {}

	/**
	 * Ctor for use with a fraction, useful for initing a variable to a non-integer number at declaration.
	 * @param numerator   The Nominator of the fraction
	 * @param denominator The Denominator of the fraction
	 */
	FixedT(int numerator, int denominator) {Raw::m_data = ((int64)numerator << Raw::dec_bits) / (int64)denominator;}

	/**
	 * Cast to a double (for outputs)
	 */
	operator double() const
	{
		return (double)Raw::m_data / (1 << Raw::dec_bits);
	}

	/**
	 * Equality operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator ==(const T &comparator) const { return Raw::m_data == (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Inequality operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator !=(const T &comparator) const { return Raw::m_data != (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Greater or equal operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator >=(const T &comparator) const { return Raw::m_data >= (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Less or equal operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator <=(const T &comparator) const { return Raw::m_data <= (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Greater than operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator >(const T &comparator) const { return Raw::m_data > (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Less than operator
	 * @param comparator The non-floating point variable we want to compare against
	 */
	template <typename T>
	bool operator <(const T &comparator) const { return Raw::m_data < (Tstorage)FixedHelperT<T>::Raw(comparator, Tdec_bits); }

	/**
	 * Addition for Fixed-point data types
	 * @param value The non-floating point number to add
	 */
	template <typename T> FixedT operator +(const T &value) const
	{
		return Raw(Raw::m_data + (Tstorage)FixedHelperT<T>::Raw(value, Tdec_bits));
	}

	/**
	 * Subtraction for Fixed-point data types
	 * @param value The non-floating point number to subtract
	 */
	template <typename T> FixedT operator -(const T &value) const
	{
		return Raw(Raw::m_data - (Tstorage)FixedHelperT<T>::Raw(value, Tdec_bits));
	}

	/**
	 * A simple multiplication for non-floating data types
	 * @param value The non floating-point factor
	 */
	template <typename T> FixedT operator *(const T &value) const
	{
		return Raw( (((int64)Raw::m_data) * FixedHelperT<T>::Raw(value)) >> FixedHelperT<T>::dec_bits);
	}

	/**
	 * A simple division for non-floating data types
	 * @param value The non floating-point divisor
	 */
	template <typename T> FixedT operator /(const T &value) const
	{
		return Raw( (((int64)Raw::m_data) << FixedHelperT<T>::dec_bits) / FixedHelperT<T>::Raw(value));
	}

	/**
	 * Addition-assignment for Fixed-point data types
	 * @param value The non-floating point number to add
	 */
	template <typename T> FixedT& operator +=(const T &value)
	{
		Raw::m_data = Raw::m_data + (Tstorage)FixedHelperT<T>::Raw(value, Tdec_bits);
		return *this;
	}

	/**
	 * Subtract-assignment for Fixed-point data types
	 * @param value The non-floating point number to subtract
	 */
	template <typename T> FixedT& operator -=(const T &value)
	{
		Raw::m_data = Raw::m_data - (Tstorage)FixedHelperT<T>::Raw(value, Tdec_bits);
		return *this;
	}

	/**
	 * Multiply-assignment for Fixed-point data types
	 * @param value The non-floating point number to multiply
	 */
	template <typename T> FixedT& operator *=(const T &value)
	{
		Raw::m_data = ((int64)Raw::m_data * (Tstorage)FixedHelperT<T>::Raw(value)) >> FixedHelperT<T>::dec_bits;
		return *this;
	}

	/**
	 * Divide-assignment for Fixed-point data types
	 * @param value The non-floating point number to use as divisor
	 */
	template <typename T> FixedT& operator /=(const T &value)
	{
		Raw::m_data = (Raw::m_data << FixedHelperT<T>::dec_bits) / (Tstorage)FixedHelperT<T>::Raw(value);
		return *this;
	}

	/**
	 * Unary minus operator
	 */
	FixedT operator -() const { return Raw(-Raw::m_data); }

	/**
	 * Stream operator, used for floating point output
	 * @param os    The stream we are going to write to
	 * @param value The Fixed-point variable we want to write in the stream
	 */
	friend std::ostream& operator << (std::ostream &os, const FixedT &value) { os << (double)value.m_data / (1ULL << Raw::dec_bits); return os; }
};

/** The value of \f$\pi\f$ with ample precision for our computations */
static const FixedT<int64, 16> PI(3141592, 1000000);
/** The number of elements used in Taylor approximations */
static const int PRECISION = 5;

/**
 * Computes a integral, positive power of a FixedT. Uses an optimized algorithm
 * to keep computational requirement at bay, by decomposing the power into powers
 * of two itself.
 *
 * @param arg   The number to compute the power for
 * @param pow   the power
 * @return      arg^pow
 * @todo        Add a nice LaTeX forumla to the documentation
 * @todo        Implement negative powers
 */
template <typename Tstorage, int Tdec_bits>
FixedT<Tstorage, Tdec_bits> pow(const FixedT<Tstorage, Tdec_bits> &arg, int pow)
{
	FixedT<int64, Tdec_bits> temp = 1;
	FixedT<int64, Tdec_bits> mul = arg;
	if (pow == 0) return (FixedT<Tstorage, Tdec_bits>)1;
	for (int i = 1; pow > 0; pow >>= 1, i <<= 1, mul *= mul) if (pow & 1) temp *= mul;

	return temp;
}

/**
 * Computes the n-th root of a number using an iterative scheme
 * @param arg  the number to compute the root from
 * @param root The "n" in the n-th root
 * @return     the root
 * @todo       Abort not after a fixed number of iterations but dynamically by step size
 * @todo       Make some better first-guess.
 */
template <typename Tstorage, int Tdec_bits>
FixedT<Tstorage, Tdec_bits> nroot(const FixedT<Tstorage, Tdec_bits> &arg, int root)
{
	FixedT<int64, 20> one = 1;
	FixedT<int64, 20> guess(6, 5);
	for (int i = 0; i < 150; i++) guess = (one / root) * ( guess * (root - 1) + arg / pow(guess, root - 1));

	return guess;
}

/**
 * Computes a single element of the Taylor series for the approximation of the cosine.
 * A single element of the Taylor series computes as:
 * \f[
 *    \frac{(-1)^nx^{2n}}{(2n)!}
 * \f]
 * In order to prevent overflows, underflows and other nasty stuff, we do not compute
 * the numerator and denominator separately in the formula but build a product of
 * fractions:
 * \f[
 *    \frac{x}{2n} * \frac{x}{2n-1} * \frac{x}{2n-2} ...
 * \f]
 * @param arg The value of "x" in above formula
 * @param pow The value of "2n" in above formula
 * @return    The value of the element of the Taylor series
 */
template <typename Tstorage, int Tdec_bits>
FixedT<Tstorage, Tdec_bits> trigbase(const FixedT<Tstorage, Tdec_bits> &arg, int pow)
{
	/* Shortcut for the very first element */
	if (pow == 0) return (FixedT<Tstorage, Tdec_bits>)1;

	/* find out the (-1)^n part, element should be negative if n is non-even
	 * but "pow" already holds 2n, account for that
	 */
	bool neg = ( (pow / 2) % 2 > 0);

	FixedT<Tstorage, Tdec_bits> element = arg / pow;
	for (--pow; pow > 0; pow --) element *= arg / pow;

	return (neg) ? -element : element;
}

/**
 * Computes the cosine of an argument using a Taylor series.
 * @param arg The number to compute the cosine from, assumes arg to be in radians.
 * @return    The cosine
 * @todo      Optimize the two while loops.
 * @note      It is possible to factor out the series to optimize even further, but is it
 *            worth the effort? It could reduce the opcount a little.
 * @see       trigbase, PRECISION
 */
template <typename Tstorage, int Tdec_bits>
FixedT<Tstorage, Tdec_bits> cos(const FixedT<Tstorage, Tdec_bits> &arg)
{
	FixedT<Tstorage, Tdec_bits> l_arg = arg;
	FixedT<Tstorage, Tdec_bits> cosine = 0;

	/* We've got to adjust the argument around zero, otherwise we will need
	 * (literally) hundreds of elements
	 */
	while (l_arg < -PI) l_arg += PI * 2;
	while (l_arg >  PI) l_arg -= PI * 2;

	/* sum up the Taylor elements */
	for(int i = 0; i < PRECISION; i++) cosine += trigbase(l_arg, 2 * i);

	return cosine;
}

/**
 * Computes the sine of an argument using a Taylor series, uses the cosine
 * computation in fact, as \f$sin(x) = cos(x - \frac{\pi}{2}))\f$
 * @param arg The number to compute the cosine from, assumes arg to be in radians.
 * @return    The sine
 */
template <typename Tstorage, int Tdec_bits>
FixedT<Tstorage, Tdec_bits> sin(const FixedT<Tstorage, Tdec_bits> &arg)
{
	return cos(arg - PI / 2);
}