#include "common.h"
#include "SampleGenerator.h"


using namespace std;

HaltonSequence SphericalSampleGenerator::sHalton;
HaltonSequence PoissonDiscSampleGenerator::sHalton;
HaltonSequence GaussianSampleGenerator::sHalton;
HaltonSequence PseudoRandomGenerator::sHalton;


SampleGenerator::SampleGenerator(int numSamples, float radius):
mNumSamples(numSamples), mRadius(radius)
{}


PoissonDiscSampleGenerator::PoissonDiscSampleGenerator(int numSamples, float radius):
SampleGenerator(numSamples, radius)
{}


void PoissonDiscSampleGenerator::Generate(float *samples) const
{
	// this is a hacky poisson sampling generator which does random dart-throwing
	// until it is not able to place any dart for a number of tries
	// in this case, the required min distance is reduced
	// the solution is a possion sampling with respect to the adjusted min distance
	// better solutions have been proposed, i.e., using hierarchical sampling

	const float maxTries = 1000;
	const float f_reduction = 0.9f;

	static HaltonSequence halton;
	float r[2];

	// generates poisson distribution on disc
	float minDist = 2.0f / sqrt((float)mNumSamples);

	//cout << "minDist before= " << minDist << endl;
	Sample2 *s = (Sample2 *)samples;

	for (int i = 0; i < mNumSamples; ++ i)
	{
		int tries = 0, totalTries = 0;

		// repeat until valid sample was found
		while (1)
		{
			++ tries;
			++ totalTries;

			halton.GetNext(2, r);

			const float rx = r[0] * 2.0f - 1.0f;
			const float ry = r[1] * 2.0f - 1.0f;

			// check if in disk, else exit early
			if (rx * rx + ry * ry > 1)
				continue;

			bool sampleValid = true;

			// check poisson property
			for (int j = 0; ((j < i) && sampleValid); ++ j)
			{
				const float dist = 
					sqrt((s[j].x - rx) * (s[j].x - rx) +
					     (s[j].y - ry) * (s[j].y - ry));
			
				if (dist < minDist)
					sampleValid = false;
			}

			if (sampleValid)
			{
				s[i].x = rx;
				s[i].y = ry;
				break;
			}

			if (tries > maxTries)
			{
				minDist *= f_reduction;
				tries = 0;
			}
		}
		
		//cout << "sample: " << samples[i].x << " " << i << " " << samples[i].y << " r: " << sqrt(samples[i].x * samples[i].x  + samples[i].y * samples[i].y) << " tries: " << totalTries << endl;
	}

	//cout << "minDist after= " << minDist << endl;
}



GaussianSampleGenerator::GaussianSampleGenerator(int numSamples, float radius):
SampleGenerator(numSamples, radius)
{}


void GaussianSampleGenerator::Generate(float *samples) const
{
	float r[2];

	const float sigma = mRadius;

	Sample2 *s = (Sample2 *)samples;

	const float p0 = 1.0f / (sigma * sqrt(2.0f * M_PI));
	const float p1 = 1.0f / (sigma * sqrt(2.0f * M_PI));

	for (int i = 0; i < mNumSamples; ++ i)
	{
		sHalton.GetNext(2, r);

		float exp_x = -(r[0] * r[0]) / (2.0f * sigma * sigma);
		float exp_y = -(r[1] * r[1]) / (2.0f * sigma * sigma);

		s[i].x = p0 * pow(M_E, exp_x);
		s[i].y = p1 * pow(M_E, exp_y);
	}

	//cout << "minDist after= " << minDist << endl;
}


PseudoRandomGenerator::PseudoRandomGenerator(int numSamples, float radius):
SampleGenerator(numSamples, radius)
{}


void PseudoRandomGenerator::Generate(float *samples) const
{
	sHalton.GetNext(2 * mNumSamples, samples);
}


SphericalSampleGenerator::SphericalSampleGenerator(int numSamples, float radius):
SampleGenerator(numSamples, radius)
{}


void SphericalSampleGenerator::Generate(float *samples) const
{
	float r[2];
	Sample3 *s = (Sample3 *)samples;
/*	
	// idea: use poisson distribution to generate samples

	PoissonDiscSampleGenerator poisson(mNumSamples, 1.0f);
	Sample2 *pSamples = new Sample2[mNumSamples];
	poisson.Generate((float *)pSamples);
*/
	for (int i = 0; i < mNumSamples; ++ i)
	{
		sHalton.GetNext(2, r);
		//r[0] = pSamples[i].x; r[1] = pSamples[i].y;

		// create stratified samples over sphere
		const float theta = 2.0f * acos(sqrt(1.0f - r[0]));
		const float phi = 2.0f * M_PI * r[1];
 
		s[i].x = mRadius * sin(theta) * cos(phi);
		s[i].y = mRadius * sin(theta) * sin(phi);
		s[i].z = mRadius * cos(theta);
	}

	//delete [] pSamples;
}