float EdgeDistance(float2 const & grad, float val)
{
	float df;
	if ((0 == grad.x()) || (0 == grad.y()))
	{
		df = 0.5f - val;
	}
	else
	{
		float2 n_grad = MathLib::abs(MathLib::normalize(grad));
		if (n_grad.x() < n_grad.y())
		{
			std::swap(n_grad.x(), n_grad.y());
		}

		float v1 = 0.5f * n_grad.y() / n_grad.x();
		if (val < v1)
		{
			df = 0.5f * (n_grad.x() + n_grad.y()) - MathLib::sqrt(2 * n_grad.x() * n_grad.y() * val);
		}
		else if (val < 1 - v1)
		{
			df = (0.5f - val) * n_grad.x();
		}
		else
		{
			df = -0.5f * (n_grad.x() + n_grad.y()) + MathLib::sqrt(2 * n_grad.x() * n_grad.y() * (1 - val));
		}
	}
	return df;
}

void IGLUShaderVariable::operator= ( float2 val )
{
	// Check for a valid shader index
	if ( m_varIdx < 0 ) return;

	// Check for type mismatches
	if ( m_isAttribute )
	{
		AssignmentToAttribute( "vec2" );
		return;
	}
	if ( m_varType != GL_FLOAT_VEC2 && m_varType != GL_DOUBLE_VEC2 )
		TypeMismatch( "vec2" );
	
	// Ensure this program is currently bound, or setting shader values fails!
	m_parent->PushProgram();

	// For types of variable that can be assigned from our input value, assign them here
	if ( m_varType == GL_FLOAT_VEC2 )
		glUniform2fv( m_varIdx, 1, val.GetConstDataPtr() );
	if ( m_varType == GL_DOUBLE_VEC2 )
		glUniform2d( m_varIdx, val.X(), val.Y() );

	// We have a short "program stack" so make sure to pop off.
	m_parent->PopProgram();
}

void InitPhysics(const float2 &gv)
{
	// FIXME: check that shaders are initialized, since renderables depend on that
	CleanPhysics();
	world = new b2World(b2Vec2(gv.x(), gv.y()));
	fixtures = new IntResourceManagerCompact<b2Fixture>([](b2Fixture *fixture)
	{
		delete (Renderable *)fixture->GetUserData();
	});
}

static inline float PolyShapeValueOnAxis(SimBody *poly, const float2 n, const float d)
{
	vector<float2> &verts = poly->transformedVertices;
	float minft = n.dot(verts[0]);

	for(u32 i=1;i<verts.size();++i)
	{
		minft = min(minft, n.dot(verts[i]));
	}

	return minft - d;
};

float PerspectiveCamera::GenerateRay( const float2& rasterSample, const float2& lensSample, Ray* ray )
{
	float3 ptRaster(rasterSample.X(), rasterSample.Y(), 0.0f);
	float3 ptCamera = TransformCoord(ptRaster, mRasterToCamera);

	ray->tMin = 0.0f;
	ray->tMax = Mathf::INFINITY;

	ray->Origin = Transform(float3(0,0,0), mCameraToWorld);
	ray->Direction = Normalize(TransformDirection(ptCamera, mCameraToWorld));
	
	return 1.0f;
}

static int FindVertsFallback(Arbiter &output_arb, SimBody *poly1, SimBody *poly2, const float2 n, const float dist)
{
	int num = 0;
	Arbiter &arb = output_arb;

	for(unsigned int i=0; i<poly1->vertices.size(); i++)
	{
		float2 v = poly1->transformedVertices[i];
		if(PolyShapeContainsVertPartial(poly2, v, n.negate()))
		{
			arb.AddContact(InitContactPoint(v, n, dist, HASH_PAIR(poly1->hashid, i)));
		}
	}

	for(unsigned int i=0; i<poly2->vertices.size(); i++)
	{
		float2 v = poly2->transformedVertices[i];
		if(PolyShapeContainsVertPartial(poly1, v, n))
		{
			arb.AddContact(InitContactPoint(v, n, dist, HASH_PAIR(poly2->hashid, i)));
		}
	}

	num = arb.numContacts;
	
	return num;
}

// set uniform to 2D vector
void Shader::SetUniform(const c8 * const name, const float2 &val)
{
	PUSH_ACTIVE_SHADER(t);
	Activate();
	glUniform2fv(GetUniformLocation(name),1, val.GetVec());
	POP_ACTIVE_SHADER(t);
};

void gui_surface::resize(const float2& buffer_size_) {
	uint2 buffer_size_abs_ = ((flags & SURFACE_FLAGS::ABSOLUTE_SIZE) == SURFACE_FLAGS::ABSOLUTE_SIZE ?
							  buffer_size_.rounded() :
							  buffer_size_ * float2(oclraster::get_width(), oclraster::get_height()));
	if(buffer.get_attachment_count() != 0 &&
	   buffer_size_abs.x == buffer_size_abs_.x && buffer_size_abs.y == buffer_size_abs_.y) {
		// same size, nothing to do here
		return;
	}
	buffer_size = buffer_size_;
	buffer_size_abs = buffer_size_abs_;
	
	delete_buffer();
	
	const bool has_depth = ((flags & SURFACE_FLAGS::NO_DEPTH) != SURFACE_FLAGS::NO_DEPTH);
	buffer = framebuffer::create_with_images(buffer_size_abs.x, buffer_size_abs.y,
											 { { IMAGE_TYPE::UINT_8, IMAGE_CHANNEL::RGBA } },
											 {
												 has_depth ? IMAGE_TYPE::FLOAT_32 : IMAGE_TYPE::NONE,
												 has_depth ? IMAGE_CHANNEL::R : IMAGE_CHANNEL::NONE
											 });
	
	// set blit vbo rectangle data
	set_offset(offset);
	
	//
	redraw();
}

float3 unproject(const float2& screen_position, float VSdepth, const float2& viewport_position, const float2& viewport_dimensions, const float4x4& proj)
{
	float2 screen_texcoord = (screen_position.xy() - viewport_position) / (viewport_dimensions - 1.0f);
	screen_texcoord.y      = 1.0f - screen_texcoord.y;
	float4 clip            = float4(VSdepth * (screen_texcoord * 2.0f - 1.0f), VSdepth * proj[2].z, VSdepth);
	float4 result          = mul(invert(proj), clip);
	return result.xyz();
}

// 2D Multi-octave Simplex noise.
//
// For each octave, a higher frequency/lower amplitude function will be added to the original.
// The higher the persistence [0-1], the more of each succeeding octave will be added.
float simplexNoise( const int octaves, const float persistence, const float scale, const float2 &v ) {
    float total = 0;
    float frequency = scale;
    float amplitude = 1;

    // We have to keep track of the largest possible amplitude,
    // because each octave adds more, and we need a value in [-1, 1].
    float maxAmplitude = 0;

    for( int i=0; i < octaves; i++ ) {
        total += simplexRawNoise( v.x() * frequency, v.y() * frequency ) * amplitude;

        frequency *= 2;
        maxAmplitude += amplitude;
        amplitude *= persistence;
    }

    return total / maxAmplitude;
}

bool SATCollide(SimBody *body1, SimBody *body2, float2 &N, f32 &t)
{
	SimBody &a = *body1;
	SimBody &b = *body2;

	if(a.vertices.size() < 2 && b.vertices.size() < 2) return false;

	Mat22 OA = Mat22::RotationMatrix(a.rotation_in_rads);
	Mat22 OB = Mat22::RotationMatrix(b.rotation_in_rads);
	Mat22 OB_T = OB.Transpose();

	Mat22 xOrient = OA * OB_T;
	float2 xOffset = (a.position - b.position) * OB_T;

	const u32 MAX_SEPERATING_AXIS = 16; 
	float2 xAxis[MAX_SEPERATING_AXIS];
	f32 taxis[MAX_SEPERATING_AXIS];
	u32 axisCount = 0;

	for(u32 i=0;i<a.seperatingAxis.size();++i)
	{
		xAxis[axisCount] = a.seperatingAxis[i] * xOrient;
		if(!IntervalIntersect(a.vertices, b.vertices, xAxis[axisCount],
			xOffset, xOrient, taxis[axisCount], t))
		{
			return false;
		}
		++axisCount;
	};
	for(u32 i=0;i<b.seperatingAxis.size();++i)
	{
		xAxis[axisCount] = b.seperatingAxis[i];
		if(!IntervalIntersect(a.vertices, b.vertices, xAxis[axisCount],
			xOffset, xOrient, taxis[axisCount], t))
		{
			return false;
		}
		++axisCount;
	};

	if(!GetMinimumTranslationVector(xAxis, taxis, axisCount, N, t))
	{
		return false;
	}

	f32 D = N.dot(xOffset);
	N =  D < 0.0f ? -N : N;

	N = N * OB;

	return true;
};

// Helper function to compute the minimal circle that contains the given three points.
// To avoid extra Sqrt() operations in hot inner loops, this function returns a Circle2D structure that has its radius squared (caller needs to compute circle.r = Sqrt(circle.r) to use)
// This is essentially a fast version of Circle2D::OptimalEnclosingCircle(a, b, c)
static Circle2D MakeCircleSq(float AB, float AC, const float2 &ab, const float2 &ac)
{
	const float AB_AC = Dot(ab,ac);
	float denom = AB*AC - AB_AC*AB_AC;
	if (Abs(denom) < 1e-5f) // Each of a, b and c lie on a straight line?
	{
		if (AB_AC > 0.f) return AB > AC ? Circle2D(ab*0.5f, AB*0.25f) : Circle2D(ac*0.5f, AC*0.25f);
		else return Circle2D((ab+ac)*0.5f, (AB + AC - 2.f*AB_AC)*0.25f);
	}
	denom = 0.5f / denom;
	float s = (AC * AB - AB_AC * AC) * denom;
	if (s < 0.f) return Circle2D(ac * 0.5f, AC * 0.25f);
	else
	{
		float t = (AC * AB - AB_AC * AB) * denom;
		if (t < 0.f) return Circle2D(ab * 0.5f, AB * 0.25f);
		else if (s + t > 1.f) return Circle2D((ab + ac) * 0.5f, (AB + AC - 2.f*AB_AC)*0.25f);
		else
		{
			const float2 center = s * ab + t * ac;
			return Circle2D(center, center.LengthSq());
		}
	}
}

bool IntervalIntersect(const std::vector<float2> &aVertices,
		const std::vector<float2> &bVertices, const float2 &axis, const float2 &relPos,
		const Mat22 &xOrient, f32 &taxis, f32 tmax)
{
	SATProjection proj0 = GetInterval(aVertices, axis * xOrient.Transpose());
	SATProjection proj1 = GetInterval(bVertices, axis);
	
	f32 h = relPos.dot(axis);
	proj0.min += h;
	proj0.max += h;

	f32 d0 = proj0.min - proj1.max;
	f32 d1 = proj1.min - proj0.max;

	if(d0 > 0.0f || d1 > 0.0f)
	{
		return false;
	}

	taxis = d0 > d1 ? d0 : d1;
	return true;
};

static inline float2 mulPerElem(const float2 &v, float f)
{
	return float2(v.getX()*f, v.getY()*f);
}

static inline float2 sumPerElem(const float2 &v, const float2& w)
{
	return float2(v.getX() + w.getX(), v.getY() + w.getY());
}

static inline uint packUnorm2x16(const float2& v)
{
	uint x = (uint)round(clamp(v.getX(), 0, 1) * 65535.0f);
	uint y = (uint)round(clamp(v.getY(), 0, 1) * 65535.0f);
	return ((uint)0x0000FFFF & x) | ((y << 16) & (uint)0xFFFF0000);
}

inline const float operator *(const float2& a, const float2& b)
{
  return a.dot(b);
}

float2 float2::Reflect(const float2 &normal) const
{
    assume(normal.IsNormalized());
    return 2.f * this->ProjectToNorm(normal) - *this;
}

void float2::Orthonormalize(float2 &a, float2 &b)
{
	assume(!a.IsZero());
	a.Normalize();
	b -= a.Dot(b) * a;
}

float float2::AngleBetweenNorm(const float2 &other) const
{
	assume(this->IsNormalized());
	assume(other.IsNormalized());
	return acos(Dot(other));
}