#ifndef Foundation_Math_misc_INCLUDED
#define Foundation_Math_misc_INCLUDED

#include <foundation/TnTypes.h>

#include <cmath>
#include <math.h>
#include <assert.h>
#include <stdlib.h>
#include <climits>

#include <foundation/Math/pe.hpp>

namespace Tn {

namespace math {

/*! For portability (i.e. because Microsquish's definitions of min,
 *  max, and abs are macros that break the Standard Template Library),
 *  min is tnmin, max is tnmax, etc. */

template <typename T>
inline T tnmin(T a, T b) { return (a<b) ? a : b; }

template <typename T>
inline T tnmax(T a, T b) { return (a>b) ? a : b; }
    
template <typename T>
inline T tnclamp(T a, T min, T max)
{
    return (a<min)?min:((a>max)?max:a);
}

template <typename T>
inline T tnabs(T a) { return (a>0) ? a : -a; }


inline bool isPowerOfTwo(unsigned int i) {
    return (i & (i - 1)) == 0;
}

template <typename T>
inline T roundToPowerOfTwo(T k)
{
    --k;

    for (size_t i = 1; i < sizeof(T) * CHAR_BIT; i <<= 1)
    {
        k = k | k >> i;
    }
    return k + 1;
}


template <typename T>
inline T degToRad(T d) {return d * T(3.14159265358979323846/180.0);}

template <typename T>
inline T radToDeg(T r) {return r * T(180.0/3.14159265358979323846);}


/*! The divmod function differs from 'div' in that it behaves
 *  correctly when either the numerator or denominator is negative.
 */

inline div_t divmod(int n, int d) {

    div_t qr = div(n, d);

    if ((qr.rem<0 && d>0) ||
        (qr.rem>0 && d<0)) {
        --qr.quot;
        qr.rem += d;
    }

    return qr;
}



/*! The mod functon differs from '%' in that it behaves correctly when
 *  either the numerator or denominator is negative.
 */

inline int mod(int n, int d) {
    int m = n % d; 
    return (((m<0) && (d>0))||
            ((m>0) && (d<0)))?(m+d):m;
}


template <typename T>
inline T wrap(T a, T min, T max);


template <>
inline float wrap(float a, float min, float max) {

    assert(max>min);

    float d = max-min;
    float s = a - min;
    float q = s/d;
    float m = q - floorf(q);
    return m * d + min;
}


template <>
inline double wrap(double a, double min, double max) {

    assert(max>min);

    double d = max-min;
    double s = a - min;
    double q = s/d;
    double m = q - floor(q);
    return m * d + min;
}


template <>
inline int32_t wrap(int32_t a, int32_t min, int32_t max) {

    int32_t d = max-min;
    int32_t s = a - min;
    int32_t r = mod(s,d);
    return r + min;
}


template <typename T>
inline T clamp(T a, T min, T max) {

    return (a<min)?min:((a>max)?max:a);
}


template <typename T>
inline T gaussian(T x, T s) {
    T c = s * sqrt(T(2)*T(3.14159265358979323846));
    return exp(-(x*x)/(T(2)*s*s))/c;
}


template <typename T>
inline T ease(T t) {
    T t_2 = t * t;
    T t_3 = t_2 * t;
    return T(3)*t_2 - T(2)*t_3;
}


template <typename U>
inline U center(const U& a, const U& b) {

    typedef typename U::component_type T;
    const int I = U::dimension;

    U result;

    pe<I-1,T>::add(&result[0], &a[0], &b[0]);
    T s = T(1)/T(2);

    pe<I-1,T>::scale(&result[0], &result[0], s);
    return result;
}


template <typename U>
inline U center(const U& a, const U& b, const U& c) {

    typedef typename U::component_type T;
    const int I = U::dimension;

    U result;

    T s = T(1)/T(3);

    pe<I-1,T>::add(&result[0], &a[0], &b[0], &c[0]);
    pe<I-1,T>::scale(&result[0], &result[0], s);


    return result;
}


template <typename U>
inline U bary(const U& a, typename U::component_type fa,
	      const U& b, typename U::component_type fb,
	      const U& c, typename U::component_type fc) {

    typedef typename U::component_type T;
    const int I = U::dimension;

    T tmp1[I];
    T tmp2[I];
    T tmp3[I];

    pe<I-1,T>::scale(tmp1, &a[0], fa);
    pe<I-1,T>::scale(tmp2, &b[0], fb);
    pe<I-1,T>::scale(tmp3, &c[0], fc);

    U result;
    pe<I-1,T>::add(&result[0], tmp1, tmp2, tmp3);
    return result;
}


template <typename T>
typename T::const_iterator nextWrap(const T& c, typename T::const_iterator i) {
    return (++i==c.end()) ? c.begin() : i;
}


template <typename T>
typename T::const_iterator prevWrap(const T& c, typename T::const_iterator i) {
    return (i==c.begin()) ? --c.end() : --i;
}


template <typename T>
typename T::iterator nextWrap(T& c, typename T::iterator i) {
    return (++i==c.end()) ? c.begin() : i;
}


template <typename T>
typename T::iterator prevWrap(T& c, typename T::iterator i) {
    return (i==c.begin()) ? --c.end() : --i;
}


template <typename T>
T recurse(T (*f)(const T&), const T& a, int i) {

    if (i<=0) return a;
    return recurse(f, f(a), i-1);
}


template <typename T>
inline T lerp(T a, T b, double s)
{
    return a + ((b - a) * s);
}


template <typename T>
inline T diffDeg(T a, T b) {

    a = wrap(a, T(-180), T(180));
    b = wrap(b, T(-180), T(180));

    return
        (a<b)
        ?
        (((b-a) < T(180)) ? (b-a) : (b - T(360) - a))
        :
        (((a-b) < T(180)) ? (b-a) : (b + T(360) - a));
}


/*!  Linearly interprets the angle via the shortest route around the
 *   circle.
 */


template <typename T>
inline T lerpDeg(T a, T b, double s) {

    a = wrap(a, T(-180), T(180));
    b = wrap(b, T(-180), T(180));

    return
        (a<b)
        ?
        (((b-a) < T(180)) ? lerp(a, b, s) : wrap(lerp(a, b - T(360), s), T(-180), T(180)))
        :
        (((a-b) < T(180)) ? lerp(a, b, s) : wrap(lerp(a, b + T(360), s), T(-180), T(180)));
}


// The Fixed point class has an explicit constructor to prevent
// unintentional mixing of number types. There is no way to define an
// explicit conversion operator in C++. So code parameterized on real
// number type that needs to convert to int needs to use an auxiliary
// function. The floor function isn't quite right because its return
// type is identical to its argument type.

inline int toInt(float a) {return int(a);}
inline int toInt(double a) {return int(a);}

inline float toFloat(float a) {return a;}
inline double toFloat(double a) {return float(a);}

// Calculates a random integer in the range [0,i) without using
// floats. Note that at least on the Win32 platform, the 'rand'
// function isn't (or at least wasn't) a very good random number
// generator.

inline int randomInt(int i) {

    assert(sizeof(long long)>=sizeof(int)*2);

    long long r = rand();
    int result = (int)((r * (long long)(i))/((long long)(RAND_MAX)+1));
    assert((result>=0) && (result<i));
    return result;
}


} // end namespace math

} // end namespace Tn


// These are templates in global namespace to provide the opportunity
// for specialization, either for optimization, or for defining the
// function for custom data types. (Like "Fixed".)

template <typename T>
inline T invSqrt(T t) {return T(1)/(T)sqrt(t);}

template <typename T>
inline T invert(const T& t) {return T(1)/t;}


#endif

