📚 课程目录第1章双三次插值实现图像放大

1. 双三次插值实现图像放大

任务说明

本案例实现双三次插值（Bicubic Interpolation）算法，将输入的低分辨率图像放大2倍，并与双线性插值（Bilinear Interpolation）进行性能对比。性能评估包含两部分：插值质量（以PSNR为指标）和计算耗时。输入图像使用标准测试图，代码将输出两种算法的放大结果及对应PSNR值。

算法实现

双三次插值利用目标像素周围4×4邻域的16个像素点，通过三次卷积核计算加权和。核函数使用经典的Keys公式（a=-0.5）：

W(x) = (a+2)|x|³ - (a+3)|x|² + 1,          |x| ≤ 1
       a|x|³ - 5a|x|² + 8a|x| - 4a,        1 < |x| < 2
       0,                                   |x| ≥ 2

实现采用分离卷积方式：先在水平方向对每行做4点加权，再在垂直方向对结果做4点加权，降低计算复杂度。边界处理采用像素复制（clamp）策略。

import numpy as np
import time
from PIL import Image
from skimage.metrics import peak_signal_noise_ratio as psnr

def cubic_kernel(x, a=-0.5):
    x = abs(x)
    if x <= 1:
        return (a + 2) * x**3 - (a + 3) * x**2 + 1
    elif x < 2:
        return a * x**3 - 5 * a * x**2 + 8 * a * x - 4 * a
    return 0

def bicubic_resize(img, scale=2):
    """双三次插值放大，img: HxWxC float32 (0-1)"""
    h, w, c = img.shape
    nh, nw = int(h * scale), int(w * scale)
    out = np.zeros((nh, nw, c), dtype=np.float32)
    # 目标像素映射到源坐标
    src_y = (np.arange(nh) + 0.5) / scale - 0.5
    src_x = (np.arange(nw) + 0.5) / scale - 0.5
    for j in range(nh):
        y = src_y[j]
        y0 = int(np.floor(y))
        dy = y - y0
        # 垂直方向4个源行索引（clamp）
        vy = np.clip(y0 - 1 + np.arange(4), 0, h - 1)
        wy = np.array([cubic_kernel(dy - (i - 1)) for i in range(4)])
        for i in range(nw):
            x = src_x[i]
            x0 = int(np.floor(x))
            dx = x - x0
            vx = np.clip(x0 - 1 + np.arange(4), 0, w - 1)
            wx = np.array([cubic_kernel(dx - (i - 1)) for i in range(4)])
            # 先水平插值4行
            tmp = np.zeros((4, c), dtype=np.float32)
            for k in range(4):
                tmp[k] = np.dot(wx, img[vy[k], vx])  # 每通道加权
            out[j, i] = np.dot(wy, tmp)
    return out

def bilinear_resize(img, scale=2):
    """双线性插值放大（对照实现）"""
    h, w, c = img.shape
    nh, nw = int(h * scale), int(w * scale)
    out = np.zeros((nh, nw, c), dtype=np.float32)
    src_y = (np.arange(nh) + 0.5) / scale - 0.5
    src_x = (np.arange(nw) + 0.5) / scale - 0.5
    for j in range(nh):
        y = src_y[j]
        y0 = int(np.floor(y))
        dy = y - y0
        y1 = min(y0 + 1, h - 1)
        for i in range(nw):
            x = src_x[i]
            x0 = int(np.floor(x))
            dx = x - x0
            x1 = min(x0 + 1, w - 1)
            top = img[y0, x0] * (1 - dx) + img[y0, x1] * dx
            bot = img[y1, x0] * (1 - dx) + img[y1, x1] * dx
            out[j, i] = top * (1 - dy) + bot * dy
    return out

# ---------- 主流程 ----------
# 加载图像并转为float32 (0-1)
img = Image.open('lena.png').convert('RGB')
img_arr = np.array(img).astype(np.float32) / 255.0

# 构造低分辨率输入（缩小再放大，模拟退化）
lr = np.array(Image.fromarray((img_arr * 255).astype(np.uint8)).resize(
    (img_arr.shape[1] // 2, img_arr.shape[0] // 2), Image.BICUBIC))
lr = lr.astype(np.float32) / 255.0

# 双三次放大
t0 = time.time()
bicubic_out = bicubic_resize(lr, scale=2)
t_bicubic = time.time() - t0

# 双线性放大
t0 = time.time()
bilinear_out = bilinear_resize(lr, scale=2)
t_bilinear = time.time() - t0

# 计算PSNR（与原始高分辨率图比较）
psnr_bicubic = psnr(img_arr, bicubic_out, data_range=1.0)
psnr_bilinear = psnr(img_arr, bilinear_out, data_range=1.0)

print(f'Bicubic  PSNR: {psnr_bicubic:.2f} dB, 耗时: {t_bicubic*1000:.1f} ms')
print(f'Bilinear PSNR: {psnr_bilinear:.2f} dB, 耗时: {t_bilinear*1000:.1f} ms')

# 保存结果
Image.fromarray((bicubic_out * 255).astype(np.uint8)).save('bicubic_2x.png')
Image.fromarray((bilinear_out * 255).astype(np.uint8)).save('bilinear_2x.png')

运行与验证

依赖安装：

pip install numpy pillow scikit-image

运行方式：将代码保存为bicubic_demo.py，与测试图像lena.png放在同一目录，执行python bicubic_demo.py。

验证判定标准：

程序正常输出两行PSNR与耗时信息，无报错；
生成bicubic_2x.png和bilinear_2x.png两张放大图像，尺寸均为原图的2倍；
双三次插值的PSNR值应高于双线性插值（通常高1~3 dB），表明重建质量更优；
双三次插值耗时通常为双线性的3~8倍，符合算法复杂度差异。