path: root/src/float/float.c

/* SPDX-License-Identifier: LGPL-2.1-only */

/*
 * The Rin Library – floating point module
 *
 * Copyright (C) 2015 Gediminas Jakutis
 */

#include "rin/float.h"
#include "rin/float_types.h"
#include "float_private.h"

/*
 * TODO:
 * · Doxygen-ize comments
 * · Diagnostics
 * · moar portability
 * · ??
 */

/*
 * signbit functions return the sign bit of a floating point number
 * in the position where it is found on the raw float data itself.
 * Set sign bit means it is negative
 */
uint32_t rin_signbitf(const float num)
{
	uint32_t ret;
	ret = rin_float_to_uint(num);
	return ret & 0x80000000u;
}

uint64_t rin_signbitd(const double num)
{
	uint64_t ret;
	ret = rin_double_to_ulong(num);
	return ret & 0x8000000000000000ull;
}

/*
 * hexstring functions format a string of the hexidecimal representation
 * of the raw floating point number data
 *
 * str has to point to enough memory to hold the said string
 *
 * mode hold modebits. Has to be either RIN_HEXSTRING_DEFAULT or any of
 * the available modes ORed togther.
 * Only one non-default mode is available at the moment.
 */

char *rin_float_to_hexstring(const float num, char *str, unsigned int mode)
{
	char buffer[16] = { '\0' };
	uint32_t tmp;
	size_t ccnt;

	tmp = rin_float_to_uint(num);
	ccnt = sprintf(buffer, "%x", tmp);

	if (mode & RIN_HEXSTRING_NOPREFIX) {
		memcpy(str, "00000000", 9);
		memcpy(&str[8 - ccnt], buffer, ccnt);
	} else {
		memcpy(str, "0x00000000", 11);
		memcpy(&str[10 - ccnt], buffer, ccnt);
	}

	return str;
}

char *rin_double_to_hexstring(const double num, char *str, unsigned int mode)
{
	char buffer[32] = { '\0' };
	uint64_t tmp;
	size_t ccnt;

	tmp = rin_double_to_ulong(num);
	ccnt = sprintf(buffer, "%lx", tmp);

	if (mode & RIN_HEXSTRING_NOPREFIX) {
		memcpy(str, "0000000000000000", 17);
		memcpy(&str[16 - ccnt], buffer, ccnt);
	} else {
		memcpy(str, "0x0000000000000000", 19);
		memcpy(&str[18 - ccnt], buffer, ccnt);
	}

	return str;
}

/*
 * These comparison functions use Units In Last Place (ULPs) to check if the numbers
 * are close enough to be considered equal, which addresses most shortcomings of
 * epsilon-based comparison.
 *
 * TODO: Still not implemented optimally and can produce results worse than epsilon
 * comparison in some situations.
 */

unsigned int rin_compare_float(const float a, const float b, uint32_t max_ulps)
{
	/* in case signs differ, "regular" comparison makes no sense */
	if (rin_signbitf(a) != rin_signbitf(b)) {
		/* unless it's a case of a signed zero, so check for that */
		if (a == b) {
			return 0;
		} else {
			return 1;
		}
	} else {
		return ((uint32_t) abs(rin_float_to_int(a) - rin_float_to_int(b))) > max_ulps;
	}
}

unsigned int rin_compare_vec2(const struct rin_vec2 a, const struct rin_vec2 b, uint32_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_float(a.x, b.x, max_ulps);
	ret += rin_compare_float(a.y, b.y, max_ulps);

	return ret;
}

unsigned int rin_compare_vec3(const struct rin_vec3 a, const struct rin_vec3 b, uint32_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_float(a.x, b.x, max_ulps);
	ret += rin_compare_float(a.y, b.y, max_ulps);
	ret += rin_compare_float(a.z, b.z, max_ulps);

	return ret;
}

unsigned int rin_compare_vec4(const struct rin_vec4 a, const struct rin_vec4 b, uint32_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_float(a.x, b.x, max_ulps);
	ret += rin_compare_float(a.y, b.y, max_ulps);
	ret += rin_compare_float(a.z, b.z, max_ulps);
	ret += rin_compare_float(a.w, b.w, max_ulps);

	return ret;
}

unsigned int rin_compare_double(const double a, const double b, uint64_t max_ulps)
{
	/* in case signs differ, "regular" comparison makes no sense */
	if (rin_signbitd(a) != rin_signbitd(b)) {
		/* unless it's a case of a signed zero, so check for that */
		if (a == b) {
			return 0;
		} else {
			return 1;
		}
	} else {
		return ((uint64_t) labs(rin_double_to_long(a) - rin_double_to_long(b))) > max_ulps;
	}
}

unsigned int rin_compare_vec2d(const struct rin_vec2d a, const struct rin_vec2d b, uint64_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_double(a.x, b.x, max_ulps);
	ret += rin_compare_double(a.y, b.y, max_ulps);

	return ret;
}

unsigned int rin_compare_vec3d(const struct rin_vec3d a, const struct rin_vec3d b, uint64_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_double(a.x, b.x, max_ulps);
	ret += rin_compare_double(a.y, b.y, max_ulps);
	ret += rin_compare_double(a.z, b.z, max_ulps);

	return ret;
}

unsigned int rin_compare_vec4d(const struct rin_vec4d a, const struct rin_vec4d b, uint64_t max_ulps)
{
	unsigned int ret;

	ret = rin_compare_double(a.x, b.x, max_ulps);
	ret += rin_compare_double(a.y, b.y, max_ulps);
	ret += rin_compare_double(a.z, b.z, max_ulps);
	ret += rin_compare_double(a.w, b.w, max_ulps);

	return ret;
}

/*
 * These [float_type]_to_[integer_type] functions do type punning in a
 * standards-compliant / without invoking undefined behavior, which
 * happens in most of the [sadly] widely used methods to do this.
 */

uint32_t rin_float_to_uint(const float num)
{
	uint32_t ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

int32_t rin_float_to_int(const float num)
{
	int32_t ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

float rin_uint_to_float(const uint32_t num)
{
	float ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

float rin_int_to_float(const int32_t num)
{
	float ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

uint64_t rin_double_to_ulong(const double num)
{
	uint64_t ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

int64_t rin_double_to_long(const double num)
{
	int64_t ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

double rin_ulong_to_double(const uint64_t num)
{
	double ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}

double rin_long_to_double(const int64_t num)
{
	double ret;

	memcpy(&ret, &num, sizeof(num));

	return ret;
}