/*
    Raydium - CQFD Corp.
    http://raydium.org/
    Released under both BSD license and Lesser GPL library license.
    See "license.txt" file.
*/

#ifndef DONT_INCLUDE_HEADERS
#include "index.h"
#else
#include "headers/math.h"
#include "headers/shadow.h"
#endif

GLfloat raydium_math_cos(GLfloat i)
{
return( (GLfloat)cos(i*PI/180) );
}

GLfloat raydium_math_sin(GLfloat i)
{
return( (GLfloat)sin(i*PI/180) );
}

GLfloat raydium_math_cos_inv(GLfloat i)
{
return(acos(i)*180/PI);
}

GLfloat raydium_math_sin_inv(GLfloat i)
{
return(asin(i)*180/PI);
}

void raydium_math_rotate(GLfloat *p, GLfloat rx, GLfloat ry, GLfloat rz, GLfloat *res)
{
res[0]= (p[0]*raydium_math_cos(ry)+(p[2]*raydium_math_cos(rx)+p[1]*raydium_math_sin(rx))*raydium_math_sin(ry))*raydium_math_cos(rz) + (p[1]*raydium_math_cos(rx)-p[2]*raydium_math_sin(rx))*raydium_math_sin(rz);
res[1]=-(p[0]*raydium_math_cos(ry)+(p[2]*raydium_math_cos(rx)+p[1]*raydium_math_sin(rx))*raydium_math_sin(ry))*raydium_math_sin(rz) + (p[1]*raydium_math_cos(rx)-p[2]*raydium_math_sin(rx))*raydium_math_cos(rz);
res[2]= (p[2]*raydium_math_cos(rx)+ p[1]*raydium_math_sin(rx))*raydium_math_cos(ry)-p[0]*raydium_math_sin(ry);
}

void raydium_math_normalize_vector4 (GLfloat * vector){
	GLfloat length;
	length=1.0/sqrtf(vector[0]*vector[0]+vector[1]*vector[1]+vector[2]*vector[2]);
	vector[0]*=length;
	vector[1]*=length;
	vector[2]*=length;
}

// pos: GLfloat[3], m: GLfloat[16]
void raydium_math_pos_to_matrix(GLfloat *pos, GLfloat *m)
{
m[0+4*0] = 1; m[0+4*1] = 0; m[0+4*2] = 0; m[0+4*3] = pos[0];
m[1+4*0] = 0; m[1+4*1] = 1; m[1+4*2] = 0; m[1+4*3] = pos[1];
m[2+4*0] = 0; m[2+4*1] = 0; m[2+4*2] = 1; m[2+4*3] = pos[2];
m[3+4*0] = 0; m[3+4*1] = 0; m[3+4*2] = 0; m[3+4*3] = 1;
}

void raydium_math_multiply_matrix4(GLfloat *matrix1,GLfloat *matrix2,GLfloat *result ){

	unsigned int i, j, k;
	for(i=0; i<4; i++)
	{
		for(j=0; j<4; j++)
		{
            GLfloat acc = 0;
		    for (k=0; k<4; k++)
		        acc += matrix1[i+k*4] * matrix2[k+j*4];

		    result[i+j*4] = acc;
		}
	}
}

// res: GLfloat[3]
void raydium_math_pos_get_modelview(GLfloat *res)
{
GLfloat tmp[16];
glGetFloatv(GL_MODELVIEW_MATRIX,tmp);
res[0]=tmp[12];
res[1]=tmp[13];
res[2]=tmp[14];
}

/* Unfinished !
// pos == res is safe
void raydium_trigo_carth_to_sphere(GLfloat *pos, GLfloat *res)
{
GLfloat r,G,T;

r=sqrt(pos[0]*pos[0] + pos[1]*pos[1] + pos[2]*pos[2]);
G=atan(pos[1]/pos[0]);
T=acos(pos[2]/r);

//printf("%f\n",r);
r=1;

if(pos[0]<0)
    res[0]=-r*cos(G);
else
    res[0]=+r*cos(G);

//res[1]=r*sin(G);
if(pos[0]<0)
    res[0]=-r*sin(T)*cos(G);
else
    res[0]=r*sin(T)*cos(G);
res[1]=r*sin(G)*cos(G);
res[2]=r*cos(T);

if(isnan(res[0])) res[0]=0;
if(isnan(res[1])) res[1]=0;
if(isnan(res[2])) res[2]=0;

}
*/

int raydium_math_pow2_next(int value)
{
int pows[]={0,2,4,8,16,32,64,128,256,512,1024,2048,4096,8192,16384,32768,65536};
int pows_count=17; // 16 + first (0)
int i;

if(value>65536 || value<0)
    {
    raydium_log("trigo: value is outside of limits of this ugly function :/");
    }

for(i=0;i<pows_count;i++)
    if(pows[i]>=value)
        return pows[i];


// should never hit this point
raydium_log("trigo: raydium_math_pow2_next: ?!!");
return -1;
}

float raydium_math_angle_from_projections(float px, float py)
{
float realangle;
//using the arccos we get the "base angle"
realangle = acos ( px );
//check quadrants
//if the Y projection is negative, the angle has to be adjusted
if ( py < 0 )
    realangle = ( float ) 2 * PI - realangle;

return realangle;
}

#ifndef RAYDIUM_NETWORK_ONLY
signed char raydium_math_point_unproject_3D(GLfloat x, GLfloat y, GLfloat z, float* resx, float* resy)
{
GLdouble sx,sy,sz;
GLdouble modelMatrix[16];
GLdouble projectionMatrix[16];
GLint   viewport[4];

raydium_camera_replace();
glGetDoublev(GL_MODELVIEW_MATRIX, modelMatrix);
glGetDoublev(GL_PROJECTION_MATRIX, projectionMatrix);
glGetIntegerv(GL_VIEWPORT, viewport);
gluProject(x,y,z,modelMatrix,projectionMatrix,viewport,&sx,&sy,&sz);

sx=(sx-viewport[0])/viewport[2]*100;
sy=(sy-viewport[1])/viewport[3]*100;
(*resx)=sx;
(*resy)=sy;

if(sz<=1.0f)