#include <stdlib.h>
#include <math.h>
#include <iostream>
#include "VecField.h"

using namespace std;


vecField::vecField(int nx, int ny, int scl)
{
	int i;
	
	parm1 = 0.2f + 0.3f*(float(rand()%1000)/1000.0f-0.5f);
	parm2 = 0.2f + 0.3f*float(rand()%1000)/1000.0f;
	parm3 = 0.1f + 0.3f*float(rand()%1000)/1000.0f;

	nY = ny;
	nX = nx;

	dX = new float*[ny];
	dY = new float*[ny];
	normX = new float*[ny];
	normY = new float*[ny];
	for(i=0;i<ny;++i)
	{
		dX[i] = new float[nx];
		dY[i] = new float[nx];
		normX[i] = new float[nx];
		normY[i] = new float[nx];
	}
	cerr << "Alloc done \n";
	updateParms();

	cerr << "Function Computed \n";
}

void vecField :: updateParms() // creates a new random field
{
	int i,j;
	float dx,dy,len;
	double lensq;


	parm1 = 0.2f + 0.3f*(float(rand()%1000)/1000.0f-0.5f);
	parm2 = 0.2f + 0.3f*float(rand()%1000)/1000.0f;
	parm3 = 0.1f + 0.3f*float(rand()%1000)/1000.0f;

	signx = 1.0f*(float(rand()%2) - 0.5f);
	
	signy = 1.0f*(float(rand()%2) - 0.5f);

	for(i=0;i<nY;++i)
	{
		for(j=0;j<nX;++j)
		{
			//getVec(float(j*scl),float(i*scl),dx,dy);  // note scale
			getVec(float(j),float(i),dx,dy); // for pacific no scale
			dX[i][j] = dx;
			dY[i][j] = dy;

			lensq = dx*dx+dy*dy;

			if (lensq > 0.001f)
			{
				len = (float)sqrt(lensq);
				
				normX[i][j] = dx/len;
				normY[i][j] = dy/len;
			}
			else
			{
				normX[i][j] = 0.0;
				normY[i][j] = 0.0;
			}
		}
	}

	for(i=0;i<nY;++i)
	{
		dX[i][0] = 0.0;
		dX[i][1] = 0.0;
		dX[i][nX-1] = 0.0;
		dX[i][nX-2] = 0.0;
		dY[i][0] = 0.0;
		dY[i][1] = 0.0;
		dY[i][nX-1] = 0.0;
		dY[i][nX-2] = 0.0;
	}
	for(i=0;i<nX;++i)
	{
		dX[0][i] = 0.0;
		dX[1][i] = 0.0;
		dX[nY-1][i] = 0.0;
		dX[nY-2][i] = 0.0;
		dY[0][i] = 0.0;
		dY[1][i] = 0.0;
		dY[nY-1][i] = 0.0;
		dY[nY-2][i] = 0.0;
	}

	double theta;

	theta = 0.75;

	sinT = (float)sin(3.14159*theta);
	cosT = (float)cos(3.14159*theta);


}


//*
void vecField::getVec(float x,float y, float &dx,float &dy)  
// this is the Differential Equation that defines the flow
// for the assignment substitute the Pacific data set
{
	float X,Y;
	int ii;

    X = x*0.15f + parm3*(x+y)*0.1f + parm3*20.0f -10.0f;

    Y = y*0.4f - 40.0f*parm1;

    dx = 0.5f*parm2  + sin(0.12f*Y) ;

	dy = 0.9f* sin(parm1*X)+ 2.5f*sin(0.1*X) +3.0f*sin(0.06f*X)*cos(Y*0.03f) +  2.7f*sin(0.015f*X)*(1.0f*0.5f*cos(Y/50.0f));


	dx *=0.5f; dy*=0.5f;
 	dx*=signx;
	dy*=signy;

	if(x > 330.0)  // the key region
	{
		dy = 0.0; 
		ii = int(y/20.0); // every 20 pixels
		if(ii>0) ii = ii - 1;
		dx = 0.08f*float( ii) +0.04f;// 0.05,0.15,0.25,2.1 ...
	};


}