#ifndef UTILS
#define UTILS

#define USE_ATAN2_APPROXIMATION

#define PI 3.14159265
#define INFINITY 1000.0
#define PHI (sqrt(5)*0.5 + 0.5)

float3x3 rotX(float a) {
	return float3x3 (1, 0, 0, 0, cos(a), -sin(a), 0, sin(a), cos(a));
}

float3x3 rotY(float a) {
	return float3x3 (cos(a), 0, sin(a), 0, 1, 0, -sin(a), 0, cos(a));
}

float3x3 rotZ(float a) {
	return float3x3 (cos(a), -sin(a), 0, sin(a), cos(a), 0, 0, 0, 1);
}

static uint rndSeed = 0;

void setRndSeed(uint seed) {
	rndSeed = seed;
}

uint hash(uint x) {
	x += x << 10, x ^= x >> 6, x += x << 3, x ^= x >> 11, x += x << 15;
	return x;
}

float rnd() {
	return asfloat((rndSeed = hash(rndSeed) & 0x007fffff) | 0x3f800000) - 1;
}

float srnd() {
	return asfloat((rndSeed = hash(rndSeed) & 0x007fffff) | 0x40000000) - 3;
}

float rnd(uint s) {
	return asfloat((s & 0x007fffff) | 0x3f800000) - 1;
}

float srnd(uint s) {
	return asfloat((s & 0x007fffff) | 0x40000000) - 3;
}

float expf(float x) {
	return asfloat(int(12102203 * x + 1064866805));
}

#ifdef USE_ATAN2_APPROXIMATION
float atan2f(float y, float x) {
	float c0 = PI * 0.25;
	float c1 = PI * 0.75;

	float absY = abs(y);
	float s = (x < 0 ? -1 : +1);
	float r = (-absY*s + x) / (x*s + absY);
	//float phi = -c0*r + (x < 0 ? c1 : c0);
	float phi = (0.1963*r*r - 0.9817)*r + (x < 0 ? c1 : c0);
	return phi * (y < 0 ? -1 : +1);
}
#else
#define atan2f atan2
#endif

float3x3 cofactorMatrix(float3x3 M) {
	float3x3 R = {
			M[1][1]*M[2][2] - M[1][2]*M[2][1], M[0][2]*M[2][1] - M[0][1]*M[2][2], M[0][1]*M[1][2] - M[0][2]*M[1][1],
			M[1][2]*M[2][0] - M[1][0]*M[2][2], M[0][0]*M[2][2] - M[0][2]*M[2][0], M[0][2]*M[1][0] - M[0][0]*M[1][2],
			M[1][0]*M[2][1] - M[1][1]*M[2][0], M[0][1]*M[2][0] - M[0][0]*M[2][1], M[0][0]*M[1][1] - M[0][1]*M[1][0]
	};
	return R;
}

float3x3 inverse(float3x3 M) {
	return cofactorMatrix(M) / determinant(M);
}

float mod(float x, float y) {
	return x - y * floor(x / y);
}

float2 mod(float2 x, float2 y) {
	return x - y * floor(x / y);
}

float3 mod(float3 x, float3 y) {
	return x - y * floor(x / y);
}

float2 mod(float2 x, float y) {
	return x - y * floor(x / y);
}

float3 mod(float3 x, float y) {
	return x - y * floor(x / y);
}

float thresholdedIdentity(float x, float threshold, float value) {
	if (x > threshold) return x;
	float u = x/threshold;
	return ((2*value - threshold) * u + 2*threshold - 3*value) * u * u + value;
}

float smoothstepQuintic(float edge0, float edge1, float x) {
	x = saturate((x - edge0) / (edge1 - edge0));
	return x*x*x*(x*(x*6 - 15) + 10);
}

float3 smoothstepQuintic(float3 edge0, float3 edge1, float3 x) {
	x = saturate((x - edge0) / (edge1 - edge0));
	return x*x*x*(x*(x*6 - 15) + 10);
}

float pulse(float center, float width, float x) {
	float t = abs(x - center);
	if (t > width) return 0;
	t /= width;
	return 1 - t*t*(3 - 2 * t);
}

float3 hsv(float h, float s, float v) {
	return lerp(1, saturate(abs(frac(h+float3(3, 2, 1)/3)*6 - 3) - 1), s) * v;
}

float3 hsvSmooth(float h, float s, float v) {
	return lerp(1, smoothstep(0, 1, saturate(abs(frac(h + float3(3, 2, 1)/3)*6 - 3) - 1)), s) * v;
}

float3 fusion(float x) {
	float t = saturate(x);
	return saturate(float3(sqrt(t), t*t*t, max(sin(3.1415*1.75*t), pow(t, 12.0))));
}

float3 uncharted2ToneMapping(float3 c) {
	float A = 0.15;
	float B = 0.50;
	float C = 0.10;
	float D = 0.20;
	float E = 0.02;
	float F = 0.30;
	float W = 11.2;

	c = ((c * (A * c + C * B) + D * E) / (c * (A * c + B) + D * F)) - E / F;
	float white = ((W * (A * W + C * B) + D * E) / (W * (A * W + B) + D * F)) - E / F;
	return c / white;
}

float3 uncharted2ToneMappingCompact(float3 c) {
	return c * (c*1.28713 + 0.153229) / (c*(c + 3.33333) + 0.4);
}

float4 hgTexture(float2 tc) {
	float4 orange = float4(1, 0.111, 0, 1);
		uint d[] = {
		0x818180ff,
		0x8f818f8f,
		0x8f8f8c8f,
		0xff808c8c
	};

	int2 c = floor(tc * 16);
	if (c.x < 0 || c.x > 15 || c.y < 0 || c.y > 15)
		return float4(0, 0, 0, 0);

	if ((d[c.y / 4u] >> ((c.y % 4u) * 8u + uint(abs(c.x - 7.5)) % 8u)) & 1u)
		return float4(1, 1, 1, 1);

	return orange;
}

float fn(float3 p) {
	float4 a = dot(floor(p), float3(1,57,21)) + float4(0,57,21,78);
	p = smoothstep(0,1,frac(p)); 
	//p = cos(frac(p)*acos(-1.))*(-.5) + .5;
	a=lerp(sin(cos(a*a)), sin(cos(a*a + 2*a + 1)), p.x);
	a.xy=lerp(a.xz,a.yw,p.y);
	return lerp(a.x,a.y,p.z);
}

float pn(float3 p) {
	float3 i = floor(p);
	float4 a = (i.x + i.y * 57 + i.z * 21); a.yzw += float3(57, 21, 78);
	float3 f = frac(p); f = f * f * (-2*f + 3);
	a = lerp(sin(cos(a) * a), sin(cos(1+a) * (1+a)), f.x);
	a.xy = lerp(a.xz, a.yw, f.y);
	return lerp(a.x, a.y, f.z);
}

float fpn(float3 p) {
	float r = pn(p * 0.06125) * 0.5;
	p.xy = cos(1) * p.xy + float2(-sin(1), sin(1)) * p.yx;
	r += pn(p * 0.125) * 0.25;
	p.zy = cos(1) * p.zy + float2(-sin(1), sin(1)) * p.yz;
	r += pn(p * 0.25) * 0.125;
	return r * (1 / (0.5 + 0.25 + 0.125));
}

#endif