插值算法

线性插值

已知两点

A(x_1, y_1)

和

B(x_2, y_2)

，以及插值参数

t (0 \leq t \leq 1）

，可以通过以下公式计算任意中间点

C(x_t, y_t)

的坐标。

x 坐标的插值：
中间点的 x 坐标 $x_t$ 是 $A$ 和 $B$ 的 $x$ 坐标按照比例 $t$ 插值的结果，公式为：
$x_t = (1 - t) \cdot x_1 + t \cdot x_2$
- 当 $t = 0$ 时， $x_t = x_1$ 。
- 当 $t = 1$ 时， $x_t = x_2$ 。
- 在 $0 \leq t \leq 1$ 的范围内， $x_t$ 线性分布在 $x_1$ 和 $x_2$ 之间。
y 坐标的插值：
中间点的 $y$ 坐标 $y_t$ 是 $A$ 和 $B$ 的 $y$ 坐标按照比例 $t$ 插值的结果，公式为：
$y_t = (1 - t) \cdot y_1 + t \cdot y_2$
- 当 $t = 0$ 时， $y_t = y_1$ 。
- 当 $t = 1$ 时， $y_t = y_2$ 。
- 在 $0 \leq t \leq 1$ 的范围内， $y_t$ 线性分布在 $y_1$ 和 $y_2$ 之间。
y 坐标的插值：
中间点的 y 坐标 $y_t$ 是 $A$ 和 $B$ 的 y 坐标按照比例 $t$ 插值的结果，公式为：
$y_t = (1 - t) \cdot y_1 + t \cdot y_2$
综合表示：
中间点 $C(x_t, y_t)$ 的坐标公式为：
$C(x_t, y_t) = \big((1 - t) \cdot x_1 + t \cdot x_2, \, (1 - t) \cdot y_1 + t \cdot y_2 \big)$

双线性插值 bilinear

双线性插值

双线性插值公式推导

双线性插值是一种常见的插值方法，主要用于二维网格中根据已知点的值来估算某个未知点的值。

1. 问题描述

给定四个已知网格点 $Q_0, Q_1, Q_2, Q_3$ 的值 $z_{00}, z_{10}, z_{01}, z_{11}$ ，以及目标插值点 $P(tx, ty)$ 的相对位置，要求插值点 $P(tx, ty)$ 的值 $z$ 。在这里刚好 $m=tx, n=ty$ ，一般来说 $tx,ty$ 是 $m,n$ 的小数部分

四个网格点的坐标分别为：
$Q_0: (0, 0), \quad Q_1: (1, 0), \quad Q_2: (0, 1), \quad Q_3: (1, 1)$
插值点 $P(tx, ty)$ 的坐标为：
$P(tx, ty), \quad 其中 tx, ty \in [0, 1]$

目标是根据 $tx, ty$ 对 $z$ 进行估算。

2. 计算 $P_1$ 和 $P_2$

步骤 1：对 $x$ 方向插值

在 $Q_0$ 和 $Q_1$ 之间插值，计算 $P_1(tx, 0)$ 的值：
$P_1(tx, 0) = z_{00} \cdot (1 - tx) + z_{10} \cdot tx$

同理，在 $Q_2$ 和 $Q_3$ 之间插值，计算 $P_2(tx, 1)$ 的值：
$P_2(tx, 1) = z_{01} \cdot (1 - tx) + z_{11} \cdot tx$

3. 根据 $P_1$ 和 $P_2$ 计算 $P(tx, ty)$

步骤 2：对 $y$ 方向插值

在 $P_1(tx, 0)$ 和 $P_2(tx, 1)$ 之间插值，计算 $P(tx, ty)$ 的值：
$P(tx, ty) = P_1(tx, 0) \cdot (1 - ty) + P_2(tx, 1) \cdot ty$

将 $P_1(tx, 0)$ 和 $P_2(tx, 1)$ 的表达式代入上式，得到：
$P(tx, ty) = \big[ z_{00} \cdot (1 - tx) + z_{10} \cdot tx \big] \cdot (1 - ty) + \big[ z_{01} \cdot (1 - tx) + z_{11} \cdot tx \big] \cdot ty$

展开化简：
$P(tx, ty) = z_{00} \cdot (1 - tx) \cdot (1 - ty) + z_{10} \cdot tx \cdot (1 - ty) + z_{01} \cdot (1 - tx) \cdot ty + z_{11} \cdot tx \cdot ty$

4. 双线性插值公式

最终公式为：
$P(tx, ty) = z_{00} \cdot (1 - tx) \cdot (1 - ty) + z_{10} \cdot tx \cdot (1 - ty) + z_{01} \cdot (1 - tx) \cdot ty + z_{11} \cdot tx \cdot ty$

5. 总结

双线性插值公式分两步：

$x$ 方向插值，计算 $P_1$ 和 $P_2$ 的值。
$y$ 方向插值，根据 $P_1$ 和 $P_2$ 的值计算目标点 $P(tx, ty)$ 。

这是一种对二维网格点的插值方法，假设 $x$ 和 $y$ 方向是线性变化的，因此公式基于线性插值的原则推导而来。

用面积理解双线性插值示意图

面积理解双线性插值

双三次插值 bicubic

双三次插值

步骤 1：计算 $x$ 方向的中间插值点 $P_0, P_1, P_2, P_3$

确定插值点和网格点
- 插值点为 $P(n, m)$ ，位于 $(n, m)$ 。
- 最近的左上角网格点为 $Q_{11} = (\lfloor n \rfloor, \lfloor m \rfloor)$ 。
- 选择 $P$ 点周围的 $4 \times 4$ 网格点：
  $q_{i,j} = (\lfloor n \rfloor + j, \lfloor m \rfloor + i), \quad i, j \in \{-1, 0, 1, 2\}$
计算 $x$ 方向的偏移
- 计算 $t_x = n - \lfloor n \rfloor$ ，即插值点在 $x$ 方向上相对于网格点的偏移。
计算 $x$ 方向的中间值 $P_0, P_1, P_2, P_3$
- 对于每一行 $i \in \{-1, 0, 1, 2\}$ ，使用 $x$ 方向的插值权重 $w_x(j, t_x)$ 计算中间插值点：
  $P_i = \sum_{j=-1}^{2} w_x(j, t_x) \cdot q_{i,j}$
- $w_x(j, t_x)$ 的权重函数 $w(x)$ 定义如下（见权重公式部分）。

步骤 2：计算 $y$ 方向的最终插值点 $P$

计算 $y$ 方向的偏移
- 计算 $t_y = m - \lfloor m \rfloor$ ，即插值点在 $y$ 方向上相对于网格点的偏移。
计算 $y$ 方向的插值值 $P$
- 使用 $y$ 方向的插值权重 $w_y(i, t_y)$ 和中间插值点 $P_0, P_1, P_2, P_3$ 计算最终插值点的值：
  $P(n, m) = \sum_{i=-1}^{2} w_y(i, t_y) \cdot P_i$

双三次插值的权重函数

权重函数 $w(x)$ 的定义为：
$w(x) = \begin{cases} 1.5 |x|^3 - 2.5 |x|^2 + 1, & \text{if } 0 \leq |x| < 1 \\ -0.5 |x|^3 + 2.5 |x|^2 - 4 |x| + 2, & \text{if } 1 \leq |x| < 2 \\ 0, & \text{if } |x| \geq 2 \end{cases}$

权重函数中的变量说明

$x$ 表示插值点与网格点之间的相对距离，定义为：
$x = j - t_x \quad \text{或} \quad x = i - t_y$
根据 $x$ 的绝对值范围，选择对应的公式计算权重。

总结步骤

计算 $x$ 方向的中间插值点
- 使用 $t_x$ 和 $w(x)$ 对每一行进行插值，得到 $P_0, P_1, P_2, P_3$ 。
计算 $y$ 方向的最终插值点
- 使用 $t_y$ 和 $w(x)$ 对 $P_0, P_1, P_2, P_3$ 进行插值，得到插值点 $P$ 。
核心公式
- $x$ 方向中间值：
  $P_i = \sum_{j=-1}^{2} w_x(j, t_x) \cdot q_{i,j}$
- $y$ 方向最终值：
  $P(n, m) = \sum_{i=-1}^{2} w_y(i, t_y) \cdot P_i$
权重函数：
$w(x) = \begin{cases} 1.5 |x|^3 - 2.5 |x|^2 + 1, & \text{if } 0 \leq |x| < 1 \\ -0.5 |x|^3 + 2.5 |x|^2 - 4 |x| + 2, & \text{if } 1 \leq |x| < 2 \\ 0, & \text{if } |x| \geq 2 \end{cases}$

cuda 仿射变换实现

双线性插值

static .,__device__ float bicubic_weight(float x) {
    const float a = -0.5f;  // 常见的值，通常取-0.5f
    if (x < 0) x = -x;
    if (x < 1) return (a + 2.0f) * powf(x, 3) - (a + 3.0f) * powf(x, 2) + 1.0f;
    if (x < 2) return a * powf(x, 3) - 5.0f * a * powf(x, 2) + 8.0f * a * x - 4.0f * a;
    return 0.0f;
}

static __global__ void warp_affine_bicubic_plane_kernel(uint8_t *src,
                                                        int src_line_size,
                                                        int src_width,
                                                        int src_height,
                                                        uint8_t *dst,
                                                        int dst_width,
                                                        int dst_height,
                                                        uint8_t const_value_st,
                                                        float *warp_affine_matrix_2_3)
{
    int dx = blockDim.x * blockIdx.x + threadIdx.x;
    int dy = blockDim.y * blockIdx.y + threadIdx.y;
    if (dx >= dst_width || dy >= dst_height)
        return;

    float m_x1 = warp_affine_matrix_2_3[0];
    float m_y1 = warp_affine_matrix_2_3[1];
    float m_z1 = warp_affine_matrix_2_3[2];
    float m_x2 = warp_affine_matrix_2_3[3];
    float m_y2 = warp_affine_matrix_2_3[4];
    float m_z2 = warp_affine_matrix_2_3[5];

    float src_x = m_x1 * dx + m_y1 * dy + m_z1;
    float src_y = m_x2 * dx + m_y2 * dy + m_z2;

    if (src_x <= 1 || src_x >= src_width - 2 || src_y <= 1 || src_y >= src_height - 2)
    {
        dst[(dy * dst_width + dx) * 3] = const_value_st;
        dst[(dy * dst_width + dx) * 3 + 1] = const_value_st;
        dst[(dy * dst_width + dx) * 3 + 2] = const_value_st;
    }
    else
    {
        // 直接在插值时获取像素值
        float wts_x[4], wts_y[4];
        for (int i = 0; i < 4; i++)
        {
            wts_x[i] = bicubic_weight(src_x - ((int)floorf(src_x) - 1 + i));
            wts_y[i] = bicubic_weight(src_y - ((int)floorf(src_y) - 1 + i));
        }

        // x 方向插值
        float p_x[4][3] = {{0}};
        for (int i = 0; i < 4; i++) 
        {
            float px_r = 0, px_g = 0, px_b = 0;
            for (int j = 0; j < 4; j++) 
            {
                int y = (int)floorf(src_y) - 1 + i;
                int x = (int)floorf(src_x) - 1 + j;
                // 边界检查并计算像素
                if (y >= 0 && y < src_height && x >= 0 && x < src_width)
                {
                    uint8_t* pixel = &src[y * src_line_size + x * 3];
                    px_r += pixel[0] * bicubic_weight(src_x - x);
                    px_g += pixel[1] * bicubic_weight(src_x - x);
                    px_b += pixel[2] * bicubic_weight(src_x - x);
                }
                else
                {
                    px_r += const_value_st * bicubic_weight(src_x - x);
                    px_g += const_value_st * bicubic_weight(src_x - x);
                    px_b += const_value_st * bicubic_weight(src_x - x);
                }
            }
            p_x[i][0] = px_r;
            p_x[i][1] = px_g;
            p_x[i][2] = px_b;
        }

        // y 方向插值并输出到目标图像
        dst[(dy * dst_width + dx) * 3]     = (uint8_t)min(max(wts_y[0] * p_x[0][0] + wts_y[1] * p_x[1][0] + wts_y[2] * p_x[2][0] + wts_y[3] * p_x[3][0] + 0.5f, 0.0f), 255.0f);
        dst[(dy * dst_width + dx) * 3 + 1] = (uint8_t)min(max(wts_y[0] * p_x[0][1] + wts_y[1] * p_x[1][1] + wts_y[2] * p_x[2][1] + wts_y[3] * p_x[3][1] + 0.5f, 0.0f), 255.0f);
        dst[(dy * dst_width + dx) * 3 + 2] = (uint8_t)min(max(wts_y[0] * p_x[0][2] + wts_y[1] * p_x[1][2] + wts_y[2] * p_x[2][2] + wts_y[3] * p_x[3][2] + 0.5f, 0.0f), 255.0f);
    }
}

双三次插值

static __device__ float bicubic_weight(float x) {
    const float a = -0.5f;  // 常见的值，通常取-0.5f
    if (x < 0) x = -x;
    if (x < 1) return (a + 2.0f) * powf(x, 3) - (a + 3.0f) * powf(x, 2) + 1.0f;
    if (x < 2) return a * powf(x, 3) - 5.0f * a * powf(x, 2) + 8.0f * a * x - 4.0f * a;
    return 0.0f;
}

static __global__ void warp_affine_bicubic_plane_kernel(uint8_t *src,
                                                        int src_line_size,
                                                        int src_width,
                                                        int src_height,
                                                        uint8_t *dst,
                                                        int dst_width,
                                                        int dst_height,
                                                        uint8_t const_value_st,
                                                        float *warp_affine_matrix_2_3)
{
    int dx = blockDim.x * blockIdx.x + threadIdx.x;
    int dy = blockDim.y * blockIdx.y + threadIdx.y;
    if (dx >= dst_width || dy >= dst_height)
        return;

    float m_x1 = warp_affine_matrix_2_3[0];
    float m_y1 = warp_affine_matrix_2_3[1];
    float m_z1 = warp_affine_matrix_2_3[2];
    float m_x2 = warp_affine_matrix_2_3[3];
    float m_y2 = warp_affine_matrix_2_3[4];
    float m_z2 = warp_affine_matrix_2_3[5];

    float src_x = m_x1 * dx + m_y1 * dy + m_z1;
    float src_y = m_x2 * dx + m_y2 * dy + m_z2;

    if (src_x <= 1 || src_x >= src_width - 2 || src_y <= 1 || src_y >= src_height - 2)
    {
        dst[(dy * dst_width + dx) * 3] = const_value_st;
        dst[(dy * dst_width + dx) * 3 + 1] = const_value_st;
        dst[(dy * dst_width + dx) * 3 + 2] = const_value_st;
    }
    else
    {
        // 直接在插值时获取像素值
        float wts_y[4];
        for (int i = 0; i < 4; i++)
        {
            wts_y[i] = bicubic_weight(src_y - ((int)floorf(src_y) - 1 + i));
        }

        // x 方向插值
        float p_x[4][3] = {{0}};
        for (int i = 0; i < 4; i++) 
        {
            float px_r = 0, px_g = 0, px_b = 0;
            for (int j = 0; j < 4; j++) 
            {
                int y = (int)floorf(src_y) - 1 + i;
                int x = (int)floorf(src_x) - 1 + j;
                // 边界检查并计算像素
                if (y >= 0 && y < src_height && x >= 0 && x < src_width)
                {
                    uint8_t* pixel = &src[y * src_line_size + x * 3];
                    px_r += pixel[0] * bicubic_weight(src_x - x);
                    px_g += pixel[1] * bicubic_weight(src_x - x);
                    px_b += pixel[2] * bicubic_weight(src_x - x);
                }
                else
                {
                    px_r += const_value_st * bicubic_weight(src_x - x);
                    px_g += const_value_st * bicubic_weight(src_x - x);
                    px_b += const_value_st * bicubic_weight(src_x - x);
                }
            }
            p_x[i][0] = px_r;
            p_x[i][1] = px_g;
            p_x[i][2] = px_b;
        }

        // y 方向插值并输出到目标图像
        dst[(dy * dst_width + dx) * 3]     = (uint8_t)min(max(wts_y[0] * p_x[0][0] + wts_y[1] * p_x[1][0] + wts_y[2] * p_x[2][0] + wts_y[3] * p_x[3][0] + 0.5f, 0.0f), 255.0f);
        dst[(dy * dst_width + dx) * 3 + 1] = (uint8_t)min(max(wts_y[0] * p_x[0][1] + wts_y[1] * p_x[1][1] + wts_y[2] * p_x[2][1] + wts_y[3] * p_x[3][1] + 0.5f, 0.0f), 255.0f);
        dst[(dy * dst_width + dx) * 3 + 2] = (uint8_t)min(max(wts_y[0] * p_x[0][2] + wts_y[1] * p_x[1][2] + wts_y[2] * p_x[2][2] + wts_y[3] * p_x[3][2] + 0.5f, 0.0f), 255.0f);
    }
}

最后编辑于：2025.01.08 15:13:24

人面猴
序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
沈念sama阅读 213,752评论 6赞 493
死咒
序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
沈念sama阅读 91,100评论 3赞 387
救了他两次的神仙让他今天三更去死
文/潘晓璐我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
开封第一讲书人阅读 159,244评论 0赞 349
道士缉凶录：失踪的卖姜人
文/不坏的土叔我叫张陵，是天一观的道长。经常有香客问我，道长，这世上最难降的妖魔是什么？我笑而不...
开封第一讲书人阅读 57,099评论 1赞 286
港岛之恋（遗憾婚礼）
正文为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
茶点故事阅读 66,210评论 6赞 385
恶毒庶女顶嫁案：这布局不是一般人想出来的
文/花漫我一把揭开白布。她就那样静静地躺着，像睡着了一般。火红的嫁衣衬着肌肤如雪。梳的纹丝不乱的头发上，一...
开封第一讲书人阅读 50,307评论 1赞 292
城市分裂传说
那天，我揣着相机与录音，去河边找鬼。笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。我是一名探鬼主播，决...
沈念sama阅读 39,346评论 3赞 412
双鸳鸯连环套：你想象不到人心有多黑
文/苍兰香墨我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
开封第一讲书人阅读 38,133评论 0赞 269
万荣杀人案实录
序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
沈念sama阅读 44,546评论 1赞 306
护林员之死
正文独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋以下内容为张勋视角年9月15日...
茶点故事阅读 36,849评论 2赞 328
白月光启示录
正文我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
茶点故事阅读 39,019评论 1赞 341
活死人
序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
沈念sama阅读 34,702评论 4赞 337
日本核电站爆炸内幕
正文年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
茶点故事阅读 40,331评论 3赞 319
男人毒药：我在死后第九天来索命
文/蒙蒙一、第九天我趴在偏房一处隐蔽的房顶上张望。院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
开封第一讲书人阅读 31,030评论 0赞 21
一桩弑父案，背后竟有这般阴谋
文/苍兰香墨我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。一阵脚步声响...
开封第一讲书人阅读 32,260评论 1赞 267
情欲美人皮
我被黑心中介骗来泰国打工，没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。一个月前我还...
沈念sama阅读 46,871评论 2赞 365
代替公主和亲
正文我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
茶点故事阅读 43,898评论 2赞 351

插值算法

线性插值

双线性插值 bilinear

双线性插值公式推导

1. 问题描述

2. 计算 和

3. 根据 和 计算

4. 双线性插值公式

5. 总结

用面积理解双线性插值示意图

双三次插值 bicubic

步骤 1：计算 方向的中间插值点

步骤 2：计算 方向的最终插值点

双三次插值的权重函数

权重函数中的变量说明

总结步骤

cuda 仿射变换实现

双线性插值

双三次插值

推荐阅读更多精彩内容

2. 计算 $P_1$ 和 $P_2$

3. 根据 $P_1$ 和 $P_2$ 计算 $P(tx, ty)$

步骤 1：计算 $x$ 方向的中间插值点 $P_0, P_1, P_2, P_3$

步骤 2：计算 $y$ 方向的最终插值点 $P$