#if !defined INCLUDE_UTILITY_ROTATION
#define INCLUDE_UTILITY_ROTATION

mat2 get_rotation_matrix(float angle) {
    float cosine = cos(angle);
    float sine = sin(angle);
    return mat2(cosine, -sine, sine, cosine);
}

// Rotation matrix using Rodriguez's rotation formula
// https://math.stackexchange.com/questions/3233842/rotation-matrix-from-axis-angle-representation
mat3 get_rotation_matrix(vec3 axis, float cosine, float sine) {
    vec3 mul = axis - axis * cosine;
    vec3 add = axis * sine;
    vec3 axis_sq = axis * axis;
    vec3 diagonal = axis_sq + (cosine - cosine * axis_sq);

    return mat3(
        diagonal.x,
        axis.x * mul.y - add.z,
        axis.x * mul.z + add.y,
        axis.x * mul.y + add.z,
        diagonal.y,
        axis.y * mul.z - add.x,
        axis.x * mul.z - mul.y,
        axis.y * mul.z + add.x,
        diagonal.z
    );
}

mat3 get_rotation_matrix(vec3 axis, float angle) {
    return get_rotation_matrix(axis, cos(angle), sin(angle));
}

#if 0
// Not tested!
// q0 + q1 i + q2 j + q3 k -> vec4(q1, q2, q3, q0)
// https://www.geeks3d.com/20141201/how-to-rotate-a-vertex-by-a-quaternion-in-glsl/

vec4 get_rotation_quaternion(vec3 axis, float angle) {
	angle *= 0.5;
	float sin_half_angle = sin(angle);
	float cos_half_angle = cos(angle);

	return vec4(axis * sin_half_angle, cos_half_angle);
}

vec4 get_rotation_quaternion(vec3 from, vec3 to) {
	vec3 axis = cross(from, to);
	float cos_angle = dot(from, to);

	// Find sine and cosine of half angle using half angle identities
	float sin_half_angle = sqrt(0.5 - 0.5 * cos_angle);
	float cos_half_angle = sqrt(0.5 + 0.5 * cos_angle);

	return vec4(axis * sin_half_angle, cos_half_angle);
}

vec3 rotate(vec4 q, vec3 v) {
	return v + 2.0 * cross(q.xyz, cross(q.xyz, v) + q.w * v);
}
#endif

#endif // INCLUDE_UTILITY_ROTATION