From 568506ecbef8352ea6127cb9a148add6cacb6c43 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jonas=20=C3=85dahl?= Date: Tue, 11 Oct 2022 11:00:55 +0200 Subject: [PATCH] cogl: Add half float implementation This implementation is copied from mesa. It uses x86_64 assembler on CPUs where it's supported, otherwise falls back on a software implementation (cogl-soft-float.c). It's intended to be used for packing and unpacking f16 format bitmaps. Part-of: --- cogl/cogl/cogl-half-float.c | 264 ++++++++++++++++++++++++++++++++++++ cogl/cogl/cogl-half-float.h | 121 +++++++++++++++++ cogl/cogl/meson.build | 2 + 3 files changed, 387 insertions(+) create mode 100644 cogl/cogl/cogl-half-float.c create mode 100644 cogl/cogl/cogl-half-float.h diff --git a/cogl/cogl/cogl-half-float.c b/cogl/cogl/cogl-half-float.c new file mode 100644 index 000000000..be2c48b81 --- /dev/null +++ b/cogl/cogl/cogl-half-float.c @@ -0,0 +1,264 @@ +/* + * Mesa 3-D graphics library + * + * Copyright (C) 1999-2007 Brian Paul All Rights Reserved. + * Copyright 2015 Philip Taylor + * Copyright 2018 Advanced Micro Devices, Inc. + * Copyright (C) 2018-2019 Intel Corporation + * + * Permission is hereby granted, free of charge, to any person obtaining a + * copy of this software and associated documentation files (the "Software"), + * to deal in the Software without restriction, including without limitation + * the rights to use, copy, modify, merge, publish, distribute, sublicense, + * and/or sell copies of the Software, and to permit persons to whom the + * Software is furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included + * in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS + * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL + * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR + * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, + * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR + * OTHER DEALINGS IN THE SOFTWARE. + */ + +/* This file comes from half_float.c in mesa. */ + +#include "config.h" + +#include + +#include "cogl/cogl-half-float.h" +#include "cogl/cogl-soft-float.h" + +typedef union +{ + float f; + int32_t i; + uint32_t u; +} FloatInt; + +/** + * \brief Rounds \c x to the nearest integer, with ties to the even integer, + * and returns the value as a long int. + */ +static inline long +cogl_lroundevenf (float x) +{ +#if defined(__SSE__) || \ + (defined(_M_IX86_FP) && (_M_IX86_FP >= 1)) || \ + (defined(_M_X64) && !defined(_M_ARM64EC)) +#if LONG_MAX == INT64_MAX + return _mm_cvtss_si64 (_mm_load_ss (&x)); +#elif LONG_MAX == INT32_MAX + return _mm_cvtss_si32 (_mm_load_ss (&x)); +#else +#error "Unsupported long size" +#endif +#else + return lrintf (x); +#endif +} + +/** + * Convert a 4-byte float to a 2-byte half float. + * + * Not all float32 values can be represented exactly as a float16 value. We + * round such intermediate float32 values to the nearest float16. When the + * float32 lies exactly between to float16 values, we round to the one with + * an even mantissa. + * + * This rounding behavior has several benefits: + * - It has no sign bias. + * + * - It reproduces the behavior of real hardware: opcode F32TO16 in Intel's + * GPU ISA. + * + * - By reproducing the behavior of the GPU (at least on Intel hardware), + * compile-time evaluation of constant packHalf2x16 GLSL expressions will + * result in the same value as if the expression were executed on the GPU. + */ +uint16_t +cogl_float_to_half_slow (float val) +{ + const FloatInt fi = {val}; + const int flt_m = fi.i & 0x7fffff; + const int flt_e = (fi.i >> 23) & 0xff; + const int flt_s = (fi.i >> 31) & 0x1; + int s, e, m = 0; + uint16_t result; + + /* sign bit */ + s = flt_s; + + /* handle special cases */ + if ((flt_e == 0) && (flt_m == 0)) + { + /* zero */ + /* m = 0; - already set */ + e = 0; + } + else if ((flt_e == 0) && (flt_m != 0)) + { + /* denorm -- denorm float maps to 0 half */ + /* m = 0; - already set */ + e = 0; + } + else if ((flt_e == 0xff) && (flt_m == 0)) + { + /* infinity */ + /* m = 0; - already set */ + e = 31; + } + else if ((flt_e == 0xff) && (flt_m != 0)) + { + /* Retain the top bits of a NaN to make sure that the quiet/signaling + * status stays the same. + */ + m = flt_m >> 13; + if (!m) + m = 1; + e = 31; + } + else { + /* regular number */ + const int new_exp = flt_e - 127; + if (new_exp < -14) + { + /* The float32 lies in the range (0.0, min_normal16) and is rounded + * to a nearby float16 value. The result will be either zero, subnormal, + * or normal. + */ + e = 0; + m = cogl_lroundevenf ((1 << 24) * fabsf (fi.f)); + } + else if (new_exp > 15) + { + /* map this value to infinity */ + /* m = 0; - already set */ + e = 31; + } + else { + /* The float32 lies in the range + * [min_normal16, max_normal16 + max_step16) + * and is rounded to a nearby float16 value. The result will be + * either normal or infinite. + */ + e = new_exp + 15; + m = cogl_lroundevenf (flt_m / (float) (1 << 13)); + } + } + + g_assert (0 <= m && m <= 1024); + if (m == 1024) + { + /* The float32 was rounded upwards into the range of the next exponent, + * so bump the exponent. This correctly handles the case where f32 + * should be rounded up to float16 infinity. + */ + ++e; + m = 0; + } + + result = (s << 15) | (e << 10) | m; + return result; +} + +uint16_t +cogl_float_to_float16_rtz_slow (float val) +{ + return cogl_float_to_half_rtz_slow (val); +} + +/** + * Convert a 2-byte half float to a 4-byte float. + * Based on code from: + * http://www.opengl.org/discussion_boards/ubb/Forum3/HTML/008786.html + */ +float +cogl_half_to_float_slow (uint16_t val) +{ + FloatInt infnan; + FloatInt magic; + FloatInt f32; + + infnan.u = 0x8f << 23; + infnan.f = 65536.0f; + magic.u = 0xef << 23; + + /* Exponent / Mantissa */ + f32.u = (val & 0x7fff) << 13; + + /* Adjust */ + f32.f *= magic.f; + /* XXX: The magic mul relies on denorms being available */ + + /* Inf / NaN */ + if (f32.f >= infnan.f) + f32.u |= 0xff << 23; + + /* Sign */ + f32.u |= (uint32_t)(val & 0x8000) << 16; + + return f32.f; +} + +/** + * Convert 0.0 to 0x00, 1.0 to 0xff. + * Values outside the range [0.0, 1.0] will give undefined results. + */ +uint8_t cogl_half_to_unorm8 (uint16_t val) +{ + const int m = val & 0x3ff; + const int e = (val >> 10) & 0x1f; + const int s = (val >> 15) & 0x1; + + /* v = round_to_nearest (1.mmmmmmmmmm * 2^(e-15) * 255) + * = round_to_nearest ((1.mmmmmmmmmm * 255) * 2^(e-15)) + * = round_to_nearest ((1mmmmmmmmmm * 255) * 2^(e-25)) + * = round_to_zero ((1mmmmmmmmmm * 255) * 2^(e-25) + 0.5) + * = round_to_zero (((1mmmmmmmmmm * 255) * 2^(e-24) + 1) / 2) + * + * This happens to give the correct answer for zero/subnormals too + */ + g_assert (s == 0 && val <= FP16_ONE); /* check 0 <= this <= 1 */ + /* (implies e <= 15, which means the bit-shifts below are safe) */ + + uint32_t v = ((1 << 10) | m) * 255; + v = ((v >> (24 - e)) + 1) >> 1; + return v; +} + +/** + * Takes a uint16_t, divides by 65536, converts the infinite-precision + * result to fp16 with round-to-zero. Used by the ASTC decoder. + */ +uint16_t cogl_uint16_div_64k_to_half (uint16_t v) +{ + /* Zero or subnormal. Set the mantissa to (v << 8) and return. */ + if (v < 4) + return v << 8; + + /* Count the leading 0s in the uint16_t */ + int n = __builtin_clz (v) - 16; + + /* Shift the mantissa up so bit 16 is the hidden 1 bit, + * mask it off, then shift back down to 10 bits + */ + int m = (((uint32_t)v << (n + 1)) & 0xffff ) >> 6; + + /* (0{n} 1 X{15-n}) * 2^-16 + * = 1.X * 2^(15-n-16) + * = 1.X * 2^(14-n - 15) + * which is the FP16 form with e = 14 - n + */ + int e = 14 - n; + + g_assert (e >= 1 && e <= 30); + g_assert (m >= 0 && m < 0x400); + + return (e << 10) | m; +} diff --git a/cogl/cogl/cogl-half-float.h b/cogl/cogl/cogl-half-float.h new file mode 100644 index 000000000..cee437497 --- /dev/null +++ b/cogl/cogl/cogl-half-float.h @@ -0,0 +1,121 @@ +/* + * Copyright (C) 1999-2007 Brian Paul All Rights Reserved. + * Copyright (C) 2018-2019 Intel Corporation + * + * Permission is hereby granted, free of charge, to any person obtaining a + * copy of this software and associated documentation files (the "Software"), + * to deal in the Software without restriction, including without limitation + * the rights to use, copy, modify, merge, publish, distribute, sublicense, + * and/or sell copies of the Software, and to permit persons to whom the + * Software is furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included + * in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS + * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL + * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR + * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, + * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR + * OTHER DEALINGS IN THE SOFTWARE. + */ + +/* This source file is originally from Mesa (src/util/half_float.h). */ + +#ifndef COGL_HALF_FLOAT_H +#define COGL_HALF_FLOAT_H + +#include "cogl/cogl-types.h" + +#include +#include + +#include "cogl/cogl-cpu-caps.h" +#include "cogl/cogl-soft-float.h" + +#ifdef __x86_64 +#include +#endif + +#define FP16_ONE ((uint16_t) 0x3c00) +#define FP16_ZERO ((uint16_t) 0) + +COGL_EXPORT +uint16_t cogl_float_to_half_slow (float val); + +COGL_EXPORT +float cogl_half_to_float_slow (uint16_t val); + +COGL_EXPORT +uint8_t cogl_half_to_unorm8 (uint16_t v); + +COGL_EXPORT +uint16_t cogl_uint16_div_64k_to_half (uint16_t v); + +COGL_EXPORT +uint16_t cogl_float_to_float16_rtz_slow (float val); + +static inline uint16_t +cogl_float_to_half (float val) +{ +#ifdef __x86_64 + if (cogl_cpu_has_cap (COGL_CPU_CAP_F16C)) + { + __m128 in = {val}; + __m128i out; + + /* $0 = round to nearest */ + __asm volatile ("vcvtps2ph $0, %1, %0" : "=v" (out) : "v" (in)); + return out[0]; + } +#endif + return cogl_float_to_half_slow (val); +} + +static inline float +cogl_half_to_float (uint16_t val) +{ +#ifdef __x86_64 + if (cogl_cpu_has_cap (COGL_CPU_CAP_F16C)) + { + __m128i in = {val}; + __m128 out; + + __asm volatile ("vcvtph2ps %1, %0" : "=v" (out) : "v" (in)); + return out[0]; + } +#endif + return cogl_half_to_float_slow (val); +} + +static inline uint16_t +cogl_float_to_float16_rtz (float val) +{ +#ifdef __x86_64 + if (cogl_cpu_has_cap (COGL_CPU_CAP_F16C)) + { + __m128 in = {val}; + __m128i out; + + /* $3 = round towards zero (truncate) */ + __asm volatile ("vcvtps2ph $3, %1, %0" : "=v" (out) : "v" (in)); + return out[0]; + } +#endif + return cogl_float_to_float16_rtz_slow (val); +} + +static inline uint16_t +cogl_float_to_float16_rtne (float val) +{ + return cogl_float_to_half (val); +} + +static inline gboolean +cogl_half_is_negative (uint16_t h) +{ + return !!(h & 0x8000); +} + +#endif /* COGL_HALF_FLOAT_H */ diff --git a/cogl/cogl/meson.build b/cogl/cogl/meson.build index cd4be09e2..ddfd40173 100644 --- a/cogl/cogl/meson.build +++ b/cogl/cogl/meson.build @@ -198,6 +198,8 @@ cogl_sources = [ 'cogl-glib-source.c', 'cogl-glsl-shader-boilerplate.h', 'cogl-graphene.c', + 'cogl-half-float.c', + 'cogl-half-float.h', 'cogl-i18n-private.h', 'cogl-index-buffer-private.h', 'cogl-index-buffer.c',