地形生成纹理滤波难题：三种插值算法仍存棋盘格，求平滑方案

我在地形生成流程中采用低分辨率地形纹理，先后尝试了双线性插值（Bilinear interpolation）、约束双三次插值（Constrained Bicubic Interpolation）、Catmull-Rom三种插值算法，但生成的地形始终存在明显棋盘格图案，希望找到能实现地形平滑效果的解决方法。

插值算法效果展示 #

• 双线性插值：
  
• 约束双三次插值：
  
• Catmull-Rom插值：
  

算法来源说明 #

前两种插值算法取自技术文章，Catmull-Rom算法取自ShaderToy，相关参考信息已嵌入源代码注释中。

算法实现代码 #

双线性插值 #

auto sampleData = [this](FVector2D __uv) {
    auto W_function = [](auto X, auto Y) {
        return X * Y;
    };

float X = __uv.X * float(gX - 1);
float Y = __uv.Y * float(gY - 1);
if (X < 0) X = 0;
if (Y < 0) Y = 0;
if (X > gX - 1) X = gX - 1;
if (Y > gY - 1) Y = gY - 1;
float PixelColor = 0, PixelColor1 = 0, PixelColor2 = 0, PixelColor3 = 0;
PixelColor = FormatedTopographyImageData[int(X) * (gY)+int(Y)];

if (X <= gX - 2 && Y <= gY - 2)
    PixelColor1 = FormatedTopographyImageData[(int(X) + 1) * (gY)+int(Y) + 1];
if (Y <= gY - 2)
    PixelColor2 = FormatedTopographyImageData[int(X) * (gY)+int(Y) + 1];
if (X <= gX - 2)
    PixelColor3 = FormatedTopographyImageData[(int(X) + 1) * (gY)+int(Y)];

float temp;
float FinalColor1 =
    W_function((1.0f - modf(X, &temp)), (1.0f - modf(Y, &temp))) * PixelColor +
    W_function(modf(X, &temp), (1.0f - modf(Y, &temp))) * PixelColor3 + //+X
    W_function((1.0f - modf(X, &temp)), modf(Y, &temp)) * PixelColor2 + //+Y
    W_function(modf(X, &temp), modf(Y, &temp)) * PixelColor1; //+X +Y

return FinalColor1;
};

约束双三次插值 #

auto sampleDataNotLinear = [this](FVector2D __uv) {
    auto W_function = [](auto X, auto Y) {
        return X * X * Y * Y 
            *( 9 - 6 * X - 6 * Y + 4 * X * Y);
    };

    float X = __uv.X * float(gX - 1);
    float Y = __uv.Y * float(gY - 1);
    if (X < 0) X = 0;
    if (Y < 0) Y = 0;
    if (X > gX - 1) X = gX - 1;
    if (Y > gY - 1) Y = gY - 1;
    float PixelColor = 0, PixelColor1 = 0, PixelColor2 = 0, PixelColor3 = 0;
    PixelColor = FormatedTopographyImageData[int(X) * (gY)+int(Y)];

    if (X <= gX - 2 && Y <= gY - 2)
        PixelColor1 = FormatedTopographyImageData[(int(X) + 1) * (gY)+int(Y) + 1];
    if (Y <= gY - 2)
        PixelColor2 = FormatedTopographyImageData[int(X) * (gY)+int(Y) + 1];
    if (X <= gX - 2)
        PixelColor3 = FormatedTopographyImageData[(int(X) + 1) * (gY)+int(Y)];

    float temp;
    float FinalColor1 =
        W_function((1.0f - modf(X, &temp)), (1.0f - modf(Y, &temp))) * PixelColor +
        W_function(modf(X, &temp), (1.0f - modf(Y, &temp))) * PixelColor3 + //+X
        W_function((1.0f - modf(X, &temp)), modf(Y, &temp)) * PixelColor2 + //+Y
        W_function(modf(X, &temp), modf(Y, &temp)) * PixelColor1; //+X +Y

    return FinalColor1;
};

Catmull-Rom插值 #

auto fract = [](float n) {
    float temp;
    return modf(n, &temp);
};
// note: (not) entirely stolen from https://gist.github.com/TheRealMJP/c83b8c0f46b63f3a88a5986f4fa982b1
//
// Samples a texture with Catmull-Rom filtering, using 9 texture fetches instead of 16.
// See http://vec3.ca/bicubic-filtering-in-fewer-taps/ for more details
auto SampleTextureCatmullRom = [this, &sampleData](FVector2D _uv)
{
    // We're going to sample a a 4x4 grid of texels surrounding the target UV coordinate. We'll do this by rounding
    // down the sample location to get the exact center of our "starting" texel. The starting texel will be at
    // location [1, 1] in the grid, where [0, 0] is the top left corner.
    FVector2D samplePos = FVector2D(_uv.X * float(gX), _uv.Y * float(gY));
    FVector2D texPos1 = FVector2D(floor(samplePos.X - 0.5f) + 0.5f, floor(samplePos.Y - 0.5f) + 0.5f);

    // Compute the fractional offset from our starting texel to our original sample location, which we'll
    // feed into the Catmull-Rom spline function to get our filter weights.
    FVector2D f = FVector2D(samplePos.X - texPos1.X, samplePos.Y - texPos1.Y);

    // Compute the Catmull-Rom weights using the fractional offset that we calculated earlier.
    // These equations are pre-expanded based on our knowledge of where the texels will be located,
    // which lets us avoid having to evaluate a piece-wise function.

    FVector2D w0 = f * (FVector2D(-0.5f) + f * (FVector2D(1.0f) - 0.5f * f));
    FVector2D w1 = FVector2D(1.0f) + f * f * (FVector2D(-2.5f) + 1.5f * f);
    FVector2D w2 = f * (FVector2D(0.5f) + f * (FVector2D(2.0f) - 1.5f * f));
    FVector2D w3 = f * f * (FVector2D(-0.5f) + 0.5f * f);

    // Work out weighting factors and sampling offsets that will let us use bilinear filtering to
    // simultaneously evaluate the middle 2 samples from the 4x4 grid.
    FVector2D w12 = w1 + w2;
    FVector2D offset12 = w2 / w12;

    // Compute the final UV coordinates we'll use for sampling the texture
    FVector2D texPos0 = texPos1 - FVector2D(1.0f);
    FVector2D texPos3 = texPos1 + FVector2D(2.0f);
    FVector2D texPos12 = texPos1 + offset12;

    texPos0.X /= float(gX);
    texPos0.Y /= float(gY);
    texPos3.X /= float(gX);
    texPos3.Y /= float(gY);
    texPos12.X /= float(gX);
    texPos12.Y /= float(gY);

    float result = 0.0f;
    result += sampleData(FVector2D(texPos0.X, texPos0.Y)) * w0.X * w0.Y;
    result += sampleData(FVector2D(texPos12.X, texPos0.Y)) * w12.X * w0.Y;
    result += sampleData(FVector2D(texPos3.X, texPos0.Y)) * w3.X * w0.Y;

    result += sampleData(FVector2D(texPos0.X, texPos12.Y)) * w0.X * w12.Y;
    result += sampleData(FVector2D(texPos12.X, texPos12.Y)) * w12.X * w12.Y;
    result += sampleData(FVector2D(texPos3.X, texPos12.Y)) * w3.X * w12.Y;

    result += sampleData(FVector2D(texPos0.X, texPos3.Y)) * w0.X * w3.Y;
    result += sampleData(FVector2D(texPos12.X, texPos3.Y)) * w12.X * w3.Y;
    result += sampleData(FVector2D(texPos3.X, texPos3.Y)) * w3.X * w3.Y;

    return FVector(result, result, result);
};

内容的提问来源于stack exchange，提问作者Kedriik