为什么我在 C++ 中的圆形抗锯齿算法会给出不对称的结果?
why does my circle anti aliasing algorithm in C++ give asymmetrical results?
一些背景:
所以我的计划是在 C++ 中创建一个点画算法,我基本上只是计划为圆的每个半径存储一大堆数据以写入 OpenGL 中的纹理贴图我不确定这是否是正确的做法,但我喜欢
会比计算机动态计算每个圆的半径更快,尤其是当许多圆的大小相同时,我的计划是创建一个函数,它只编写一个充满半径的整个文本文档,达到一定大小并且该数据将按位存储在长数组 std::array <long> bit = {0x21, 0x0A ect... }
中,以便我可以对 4X4 值数组进行编码,并将 2 位分配给每个像素的抗锯齿值,但是要创建这个抗锯齿圆数据库,我需要编写一个我经常出错的函数;
真题:
所以这可能看起来很懒,但我可以保证我已经尽了一切努力来解决我在这里犯的错误基本上我已经通过将像素分成子像素将这段代码编写为 anti=alias 但是看起来返回大于 1 的值,这是不可能的,因为我已将每个像素分成 100 个大小为 0.01
的像素
float CircleConst::PixelAA(int I, int J)
{
float aaValue = 0;
for (float i = (float) I; i < I + 1; i += 0.1f)
{
for (float j = (float) J; j < J + 1; j += 0.1f)
{
if ((pow((i - center), 2) + pow((j - center), 2) < pow(rad, 2)))
aaValue += 0.01f;
}
}
return aaValue;
}
这里还有写实际圆的代码
CircleConst::CircleConst(float Rad)
{
rad = Rad;
dataSize = (unsigned int) ceil(2 * rad);
center = (float) dataSize/2;
arrData.reserve((int) pow(dataSize, 2));
for (int i = 0; i < dataSize; i++)
{
for (int j = 0; j < dataSize; j++)
{
if ( CircleBounds(i, j, rad-1) )
arrData.push_back(1);
else if (!CircleBounds(i, j, rad - 1) && CircleBounds(i, j, rad + 1))
{
arrData.push_back(PixelAA(i,j));
}
else
arrData.push_back(0);
}
}
}
所以我注意到在没有抗锯齿的情况下,圆的书写方式移动了一行,但这可以通过将圆心的值更改为dataSize/2 - 0.5f
来解决,但这会导致问题稍后当圆与抗锯齿不对称时,这里是半径 3.5
的示例
0.4 1.0 1.1 1.1 1.1 0.4 0.0
1.0 1.0 1.0 1.0 1.0 1.1 0.2
1.1 1.0 1.0 1.0 1.0 1.0 0.5
1.1 1.0 1.0 1.0 1.0 1.0 0.5
1.1 1.0 1.0 1.0 1.0 1.0 0.2
0.4 1.1 1.0 1.0 1.0 0.5 0.0
0.0 0.2 0.5 0.5 0.2 0.0 0.0
如您所见,有些值超过了 1.0,这应该是不可能的,我确信对于为什么会这样有一个明显的答案,但我完全没有找到它。
问题出在这样的行上:
for (float i = (float) I; i < I + 1; i += 0.1f)
浮点数不能无限精确地存储或操作。通过将一个浮点数重复添加到另一个浮点数,误差会累积。这就是您看到值高于 1.0 的原因。
解决方案是使用整数类型进行迭代并计算所需的浮点数。例如:
for (unsigned i = 0U; i < 10U; ++i)
{
float x = 0.1F * static_cast<float>(i);
printf("%f\n", x);
}
除了@Yun(浮点数的舍入误差)表示的,还要注意采样点(必须在像素中心)
这是您的代码,经过一些修改和添加:
#include <iostream>
#include <vector>
#include <iomanip>
#include <math.h>
float rad, radSquared, center;
const int filterSize = 8;
const float invFilterSize = 1.0f / filterSize;
// Sample the circle returning 1 when inside, 0 otherwise.
int SampleCircle(int i, int j) {
float di = (i + 0.5f) * invFilterSize - center;
float dj = (j + 0.5f) * invFilterSize - center;
return ((di * di + dj * dj) < radSquared) ? 1 : 0;
}
// NOTE: This sampling method works with any filter size.
float PixelAA(int I, int J)
{
int aaValue = 0;
for (int i = 0; i < filterSize; ++i)
for (int j = 0; j < filterSize; ++j)
aaValue += SampleCircle(I + i, J + j);
return (float)aaValue / (float)(filterSize * filterSize);
}
// NOTE: This sampling method works only with filter sizes that are power of two.
float PixelAAQuadTree(int i, int j, int filterSize)
{
if (filterSize == 1)
return (float)SampleCircle(i, j);
// We sample the four corners of the filter. Note that on left and bottom corners
// 1 is subtracted to avoid sampling overlap.
int topLeft = SampleCircle(i, j);
int topRight = SampleCircle(i + filterSize - 1, j);
int bottomLeft = SampleCircle(i, j + filterSize - 1);
int bottomRight = SampleCircle(i + filterSize - 1, j + filterSize - 1);
// If all samples have same value we can stop here. All samples lies outside or inside the circle.
if (topLeft == topRight && topLeft == bottomLeft && topLeft == bottomRight)
return (float)topLeft;
// Half the filter dimension.
filterSize /= 2;
// Recurse.
return (PixelAAQuadTree(i, j, filterSize) +
PixelAAQuadTree(i + filterSize, j, filterSize) +
PixelAAQuadTree(i, j + filterSize, filterSize) +
PixelAAQuadTree(i + filterSize, j + filterSize, filterSize)) / 4.0f;
}
void CircleConst(float Rad, bool useQuadTree)
{
rad = Rad;
radSquared = rad * rad;
center = Rad;
int dataSize = (int)ceil(rad * 2);
std::vector<float> arrData;
arrData.reserve(dataSize * dataSize);
if (useQuadTree)
{
for (int i = 0; i < dataSize; i++)
for (int j = 0; j < dataSize; j++)
arrData.push_back(PixelAAQuadTree(i * filterSize, j * filterSize, filterSize));
}
else
{
for (int i = 0; i < dataSize; i++)
for (int j = 0; j < dataSize; j++)
arrData.push_back(PixelAA(i * filterSize, j * filterSize));
}
for (int i = 0; i < dataSize; i++)
{
for (int j = 0; j < dataSize; j++)
std::cout << std::fixed << std::setw(2) << std::setprecision(2)
<< std::setfill('0') << arrData[i + j * dataSize] << " ";
std::cout << std::endl;
}
}
int main() {
CircleConst(3.5f, false);
std::cout << std::endl;
CircleConst(4.0f, false);
std::cout << std::endl;
std::cout << std::endl;
CircleConst(3.5f, true);
std::cout << std::endl;
CircleConst(4.0f, true);
return 0;
}
给出这些结果(第二个使用四叉树优化计算 AA 值所需的样本数):
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.84 1.00 1.00 1.00 1.00 1.00 0.84
1.00 1.00 1.00 1.00 1.00 1.00 1.00
0.84 1.00 1.00 1.00 1.00 1.00 0.84
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.84 1.00 1.00 1.00 1.00 1.00 0.84
1.00 1.00 1.00 1.00 1.00 1.00 1.00
0.84 1.00 1.00 1.00 1.00 1.00 0.84
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
进一步说明:
- 你可以在这里看到四叉树是如何工作的https://en.wikipedia.org/wiki/Quadtree
- 您可以进一步修改代码并实现定点数学 (https://en.wikipedia.org/wiki/Fixed-point_arithmetic),它没有像浮点数那样的舍入问题,因为数字总是表示为整数
- 鉴于这些数据是预计算阶段的一部分,代码的简单性优先于性能
一些背景:
所以我的计划是在 C++ 中创建一个点画算法,我基本上只是计划为圆的每个半径存储一大堆数据以写入 OpenGL 中的纹理贴图我不确定这是否是正确的做法,但我喜欢
会比计算机动态计算每个圆的半径更快,尤其是当许多圆的大小相同时,我的计划是创建一个函数,它只编写一个充满半径的整个文本文档,达到一定大小并且该数据将按位存储在长数组 std::array <long> bit = {0x21, 0x0A ect... }
中,以便我可以对 4X4 值数组进行编码,并将 2 位分配给每个像素的抗锯齿值,但是要创建这个抗锯齿圆数据库,我需要编写一个我经常出错的函数;
真题:
所以这可能看起来很懒,但我可以保证我已经尽了一切努力来解决我在这里犯的错误基本上我已经通过将像素分成子像素将这段代码编写为 anti=alias 但是看起来返回大于 1 的值,这是不可能的,因为我已将每个像素分成 100 个大小为 0.01
的像素float CircleConst::PixelAA(int I, int J)
{
float aaValue = 0;
for (float i = (float) I; i < I + 1; i += 0.1f)
{
for (float j = (float) J; j < J + 1; j += 0.1f)
{
if ((pow((i - center), 2) + pow((j - center), 2) < pow(rad, 2)))
aaValue += 0.01f;
}
}
return aaValue;
}
这里还有写实际圆的代码
CircleConst::CircleConst(float Rad)
{
rad = Rad;
dataSize = (unsigned int) ceil(2 * rad);
center = (float) dataSize/2;
arrData.reserve((int) pow(dataSize, 2));
for (int i = 0; i < dataSize; i++)
{
for (int j = 0; j < dataSize; j++)
{
if ( CircleBounds(i, j, rad-1) )
arrData.push_back(1);
else if (!CircleBounds(i, j, rad - 1) && CircleBounds(i, j, rad + 1))
{
arrData.push_back(PixelAA(i,j));
}
else
arrData.push_back(0);
}
}
}
所以我注意到在没有抗锯齿的情况下,圆的书写方式移动了一行,但这可以通过将圆心的值更改为dataSize/2 - 0.5f
来解决,但这会导致问题稍后当圆与抗锯齿不对称时,这里是半径 3.5
0.4 1.0 1.1 1.1 1.1 0.4 0.0
1.0 1.0 1.0 1.0 1.0 1.1 0.2
1.1 1.0 1.0 1.0 1.0 1.0 0.5
1.1 1.0 1.0 1.0 1.0 1.0 0.5
1.1 1.0 1.0 1.0 1.0 1.0 0.2
0.4 1.1 1.0 1.0 1.0 0.5 0.0
0.0 0.2 0.5 0.5 0.2 0.0 0.0
如您所见,有些值超过了 1.0,这应该是不可能的,我确信对于为什么会这样有一个明显的答案,但我完全没有找到它。
问题出在这样的行上:
for (float i = (float) I; i < I + 1; i += 0.1f)
浮点数不能无限精确地存储或操作。通过将一个浮点数重复添加到另一个浮点数,误差会累积。这就是您看到值高于 1.0 的原因。
解决方案是使用整数类型进行迭代并计算所需的浮点数。例如:
for (unsigned i = 0U; i < 10U; ++i)
{
float x = 0.1F * static_cast<float>(i);
printf("%f\n", x);
}
除了@Yun(浮点数的舍入误差)表示的,还要注意采样点(必须在像素中心)
这是您的代码,经过一些修改和添加:
#include <iostream>
#include <vector>
#include <iomanip>
#include <math.h>
float rad, radSquared, center;
const int filterSize = 8;
const float invFilterSize = 1.0f / filterSize;
// Sample the circle returning 1 when inside, 0 otherwise.
int SampleCircle(int i, int j) {
float di = (i + 0.5f) * invFilterSize - center;
float dj = (j + 0.5f) * invFilterSize - center;
return ((di * di + dj * dj) < radSquared) ? 1 : 0;
}
// NOTE: This sampling method works with any filter size.
float PixelAA(int I, int J)
{
int aaValue = 0;
for (int i = 0; i < filterSize; ++i)
for (int j = 0; j < filterSize; ++j)
aaValue += SampleCircle(I + i, J + j);
return (float)aaValue / (float)(filterSize * filterSize);
}
// NOTE: This sampling method works only with filter sizes that are power of two.
float PixelAAQuadTree(int i, int j, int filterSize)
{
if (filterSize == 1)
return (float)SampleCircle(i, j);
// We sample the four corners of the filter. Note that on left and bottom corners
// 1 is subtracted to avoid sampling overlap.
int topLeft = SampleCircle(i, j);
int topRight = SampleCircle(i + filterSize - 1, j);
int bottomLeft = SampleCircle(i, j + filterSize - 1);
int bottomRight = SampleCircle(i + filterSize - 1, j + filterSize - 1);
// If all samples have same value we can stop here. All samples lies outside or inside the circle.
if (topLeft == topRight && topLeft == bottomLeft && topLeft == bottomRight)
return (float)topLeft;
// Half the filter dimension.
filterSize /= 2;
// Recurse.
return (PixelAAQuadTree(i, j, filterSize) +
PixelAAQuadTree(i + filterSize, j, filterSize) +
PixelAAQuadTree(i, j + filterSize, filterSize) +
PixelAAQuadTree(i + filterSize, j + filterSize, filterSize)) / 4.0f;
}
void CircleConst(float Rad, bool useQuadTree)
{
rad = Rad;
radSquared = rad * rad;
center = Rad;
int dataSize = (int)ceil(rad * 2);
std::vector<float> arrData;
arrData.reserve(dataSize * dataSize);
if (useQuadTree)
{
for (int i = 0; i < dataSize; i++)
for (int j = 0; j < dataSize; j++)
arrData.push_back(PixelAAQuadTree(i * filterSize, j * filterSize, filterSize));
}
else
{
for (int i = 0; i < dataSize; i++)
for (int j = 0; j < dataSize; j++)
arrData.push_back(PixelAA(i * filterSize, j * filterSize));
}
for (int i = 0; i < dataSize; i++)
{
for (int j = 0; j < dataSize; j++)
std::cout << std::fixed << std::setw(2) << std::setprecision(2)
<< std::setfill('0') << arrData[i + j * dataSize] << " ";
std::cout << std::endl;
}
}
int main() {
CircleConst(3.5f, false);
std::cout << std::endl;
CircleConst(4.0f, false);
std::cout << std::endl;
std::cout << std::endl;
CircleConst(3.5f, true);
std::cout << std::endl;
CircleConst(4.0f, true);
return 0;
}
给出这些结果(第二个使用四叉树优化计算 AA 值所需的样本数):
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.84 1.00 1.00 1.00 1.00 1.00 0.84
1.00 1.00 1.00 1.00 1.00 1.00 1.00
0.84 1.00 1.00 1.00 1.00 1.00 0.84
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.84 1.00 1.00 1.00 1.00 1.00 0.84
1.00 1.00 1.00 1.00 1.00 1.00 1.00
0.84 1.00 1.00 1.00 1.00 1.00 0.84
0.36 1.00 1.00 1.00 1.00 1.00 0.36
0.00 0.36 0.84 1.00 0.84 0.36 0.00
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.97 1.00 1.00 1.00 1.00 1.00 1.00 0.97
0.70 1.00 1.00 1.00 1.00 1.00 1.00 0.70
0.16 0.95 1.00 1.00 1.00 1.00 0.95 0.16
0.00 0.16 0.70 0.97 0.97 0.70 0.16 0.00
进一步说明:
- 你可以在这里看到四叉树是如何工作的https://en.wikipedia.org/wiki/Quadtree
- 您可以进一步修改代码并实现定点数学 (https://en.wikipedia.org/wiki/Fixed-point_arithmetic),它没有像浮点数那样的舍入问题,因为数字总是表示为整数
- 鉴于这些数据是预计算阶段的一部分,代码的简单性优先于性能