/*
 * =====================================================================================
 *
 *       Filename:  jNormal.cpp
 *
 *    Description:  implementation for jNormal
 *
 *        Version:  1.0
 *        Created:  31/05/2011 22:12:17
 *       Revision:  none
 *       Compiler:  gcc
 *
 *         Author:  John Ryland (jryland), jryland@xiaofrog.com
 *        Company:  InvertedLogic
 *
 * =====================================================================================
 */

#include <jNormal.h>


#define sqr(x)      ((x)*(x))


jNormal::jNormal()
{
    m_flags = DIRTY;
}


jNormal::jNormal(jVector3f& a_vector)
{
    m_flags = DIRTY;
    m_vec = a_vector;
}


jNormal::jNormal(jNormal& a_normal)
{
    m_flags = a_normal.m_flags;
    m_vec = a_normal.m_vec;
}


void jNormal::setNormal(jNormal& a_normal)
{
    m_flags = a_normal.m_flags;
    m_vec = a_normal.m_vec;
}


void jNormal::normalize()
{
    if (m_flags == NORMALIZED)
    {
        return;
    }
    jFloat32 one_over_dist = 1.0 / sqrt( m_vec.dotProduct(m_vec) );
    m_vec *= one_over_dist; 
    m_flags = NORMALIZED;
}


void jNormal::setVector(const jVector3f& a_vector, bool a_normalized)
{
    m_flags = (a_normalized) ? 0 : DIRTY;
    m_vec = a_vector;
}


jVector3f jNormal::getVector()
{
    return m_vec;
}


jNormalIterator::jNormalIterator(jNormal& a_startNormal, jNormal& a_endNormal, int a_steps)
{
    a_startNormal.normalize();
    a_endNormal.normalize();

    const jVector3f& A = a_startNormal.getVector();
    const jVector3f& B = a_endNormal.getVector();

    m_current = A;

    jFloat32 one_over_steps = 1.0f / a_steps;
    m_delta = (B - A) * one_over_steps;

    m_distSqrForwardDiff = 2.0f * (B.sumComponents() - A.sumComponents()) * one_over_steps;
    m_distSqr = A.dotProduct(A);
}


jNormal jNormalIterator::step(int a_times)
{
    jNormal out;

    m_distSqr += m_distSqrForwardDiff * a_times;

    jFloat32 one_over_dist = 1.0 / sqrt( m_distSqr );
    
    m_current += m_delta;

    out.setVector(m_current * one_over_dist, true);

    return out;
}


jNormalIterator::jNormalIterator(jNormal& a_startNormal, jNormal& a_endNormal, jNormal& a_dotProductWith, int a_steps)
{
    const jVector3f& A = a_startNormal.getVector();
    const jVector3f& B = a_endNormal.getVector();
    jVector3f delta = (B - A) * a_dotProductWith.getVector();
    m_unNormalizedDotForwardDiff = delta.sumComponents() / a_steps;
    m_currentDot = A.dotProduct(a_dotProductWith.getVector());
    // intergerize these values
    // TODO


}


jNormalIterator::jNormalIterator(jNormal& a_startNormal, jNormal& a_endNormal, jVector3f& a_start, jVector3f& a_end, jVector3f& a_dotProductWith, int a_steps)
{
    // TODO
    const jVector3f& A = a_startNormal.getVector();
    const jVector3f& B = a_endNormal.getVector();
    jVector3f delta = (B - A) * a_dotProductWith;;
    m_unNormalizedDotForwardDiff = delta.sumComponents() / a_steps;
   
    // other is vector from point on surface to point of light
    // point on surface changes as we iterate, so we need a point A and point B for this
    // so need 2 more arguments

    m_surfacePoint = a_start;
    m_surfacePointDelta = (a_end - a_start) / a_steps;

    //m_currentDot = A.DotProduct(other);
}


jFloat32 jNormalIterator::stepDotPointF()
{
    m_surfacePoint += m_surfacePointDelta;
    jVector3f light;
    jVector3f direction = light - m_surfacePoint;
    light.normalize();
    jVector3f light2 = light.normalized();
    jMatrix3x3f mat;
    mat.normalize();
    direction /= direction.length();
    return 0.0f;
}


// Fairly portable, and fairly accurate, reasonably fast
jFloat32 jFastInvSqrt(jFloat32 number)
{
    const float threehalfs = 1.5F;

    jFloat32 x2 = number * 0.5F;
    jFloat32 y  = number;
    jInt32   i  = *(long *)&y;
    //i  = 0x5f3759df - ( i >> 1 );
    i  = 0x5f37642f - ( i >> 1 ); // more accurate constant
    y  = *(float *)&i;
    y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
    //    y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

    return y;
}


jInt8 jNormalIterator::stepDotInfinityI()
{


    // if sqrt(x) gives values from 0->255, then x is 0->65535
    // so could use a table from 0->65535 to give us our data

/*
    0  0    8  2   16  4   24  4
    1  1    9  3   17  4   25  5
    2  1   10  3   18  4   26  5
    3  1   11  3   19  4   27  5
    4  2   12  3   20  4   29  5
    5  2   13  3   21  4   30  5
    6  2   14  3   22  4   31  5
    7  2   15  3   23  4   32  5

    first 256 will all be 4bits (0-15)


    1 3 5 7 9 ...

    25 -> 5

    next will be 25 + (5*2+1) -> 36


    sqrt(100)        10  +  0 / 20
    sqrt(99)         9   + 17 / 18
    sqrt(98)         9   + 16 / 18
    sqrt(..)         9
    sqrt(81)         9   +  0 / 18
    sqrt(80)         8   + 15 / 16
    sqrt(..)         8
    sqrt(64)         8   +  0 / 16
    sqrt(63)         7   + 13 / 14

*/

    m_distSqr += m_distSqrForwardDiff;
    jFloat32 one_over_dist = 1.0 / sqrt( m_distSqr ); // SSE instruction RSQRTSS or fastInvSqrt
    m_currentDot += m_unNormalizedDotForwardDiff;
    return m_currentDot * one_over_dist;
}


jFloat32 jNormalIterator::stepDotInfinityF()
{
    m_distSqr += m_distSqrForwardDiff;
    jFloat32 one_over_dist = 1.0 / sqrt( m_distSqr ); // SSE instruction RSQRTSS or fastInvSqrt
    m_currentDot += m_unNormalizedDotForwardDiff;
    return m_currentDot * one_over_dist;

    //      1 / sqrt(x) = sqrt(x) / x

    // dot squared
    // return m_currentDot * m_currentDot / m_distSqr;

    /*
       √(a*b) = √a * √b       if a,b ≥ 0
       ln√a = (0.5)*ln(a)     if a > 0
       √(1/a) = 1 / √a        if a > 0
       √(a^2) = |a|
       (√a)^2 = a             if a ≥ 0
       1/(√a+√b) = (√a-√b)/(a-b), if a,b ≥ 0, a≠b
       1/(√a-√b) = (√a+√b)/(a-b), if a,b ≥ 0, a≠b
       |a + i*b| = √(a2+b2)
    */

    /*
        1*1 = 1      1 = 1
        2*2 = 4      1+3 = 4
        3*3 = 9      1+3+5 = 9
        4*4 = 16     1+3+5+7 = 16
        5*5 = 25     1+3+5+7+9 = 25
    */
}