我正在做一个我做RGB转换为luma的项目,而且我对-mno-sse2标志有一些舍入问题:
这是测试代码:
#include <stdio.h>
#include <stdint.h>
static double rec709_luma_coeff[3] = {0.2126, 0.7152, 0.0722};
int main()
{
uint16_t n = 242 * rec709_luma_coeff[0] + 242 * rec709_luma_coeff[1] + 242 * rec709_luma_coeff[2];
printf("%u\n", n);
return 0;
}
而这就是我得到的:
user@gentoo>gcc -mno-sse2 test.c -o test && ./test
241
user@gentoo> gcc test.c -o test && ./test
242
我想gcc使用sse2优化double
乘法,但我没有得到的是为什么优化版本是正确的。
另外,你建议我用什么来获得更一致的结果,ceil()
或floor()
?
TL:DR使用lrint(x)
或(int)rint(x)
从float转换为int,使用round-to-nearest而不是truncation。不幸的是,并非所有编译器都能有效地内联相同的数学函数。见round() for float in C++
gcc -mno-sse2
必须使用x87用于double
,即使在64位代码中也是如此。 x87寄存器的内部精度为80位,但SSE2本身在XMM寄存器中使用IEEE binary64 (aka double
)格式,因此所有临时值都会在每一步舍入为64位double
。
问题并不像the double rounding problem那样有趣(80位 - > 64位,然后是整数)。当将临时存储器存储到内存时,它也不是来自gcc -O0
(默认:没有额外的优化)舍入,因为你在一个C语句中完成了所有操作,所以它只对整个表达式使用x87寄存器。
简单地说,80位精度导致结果刚好低于242.0并被C的float-> int语义截断为241,而SSE2产生的结果恰好高于242.0,截断为242.对于x87,舍入到下一个更低对于从1到65535的任何输入,整数一致地发生,而不仅仅是242.(我使用atoi(argv[1])
制作了一个版本的程序,所以我可以测试其他值,并使用-O3
)。
请记住,int foo = 123.99999
是123,因为C使用“截断”舍入模式(朝向零)。对于非负数,这与floor
(向-Infinity一致)相同。 https://en.wikipedia.org/wiki/Floating-point_arithmetic#Rounding_modes。
double
不能完全代表系数:我用gdb
打印它们并得到:{0.21260000000000001, 0.71519999999999995, 0.0722}
。这些十进制表示可能不是base-2浮点值的精确表示。但他们足够接近,看到系数加起来0.99999999999999996
(使用任意精度计算器)。
由于x87内部精度高于系数的精度,因此我们得到四舍五入,因此n * rec709_luma_coeff[0]
中的和舍入误差等等,并且在总结结果时,~2^11
小于系数之和的差值。 1.0。 (64位有效数与53位)。
真正的问题是SSE2版本如何成功运作!大概到最接近甚至临时的圆形在足够的情况下恰好向上移动,至少在242.它恰好产生原始输入的更多情况,但它产生输入-1为5,7,10,13, 14,20 ......(来自1..1000的前1000个数字中的252个被SSE2版本“淹没”,所以它不像它总是有效。)
使用-O3
作为源,它在编译时以扩展的精度进行计算并产生精确的结果。即它汇编与printf("%u\n", n);
相同。
顺便说一句,你应该使用static
const
为你的常数,所以gcc可以更好地优化。但是,static
比普通全局要好得多,因为编译器可以看到编译单元中没有任何内容写入值或将其地址传递到任何地方,因此它可以将它们视为const
。