数字图像处理（十三）图像放大及双线性插值算法

导言

图像放大是日常学习中经常要用到的两个算法，我们首先讨论缩放的流程以及放大时如何优化双线性插值算法。
采用国际标准测试图像Lena，为了方便，我们将读入的彩色图转为灰度图进行缩放。

图像放大

和图像缩小不同，图像放大是小数据量到大数据量的处理过程，因此需要对许多未知的数据进行估计。
如果一幅$W\times H$图像要放大$k_1\times k_2$（即行放大$k_1$倍，列放大$k_1$倍），则放大后的图像大小为$int(W\ast k_1)\times int(H\ast k_1)$，这里取int的原因是乘出来的的数值可能为小数。
具体过程如下图所示：
我们有一个图像大小为4*4，每一个像素值都是1。

它的坐标矩阵如下所示

现在我们将原图（下面称为f）的行和列都放大1.5倍($k_1=1.5$，$k_2=1.5$)。则放大后的图像（下面称为g）的 widht=6=$4\times 1.5$，height=6=$4\times 1.5$。接着我们要计算g中的坐标在f中映射。如下图所示：

我们计算出g中的(0,1)点应该对应f中的(0,0.67)点，但是f中无此坐标，所有我们使用线性插值法来计算(0,0.67)处的像素值。
线性插值法如果你没有了解过，建议先去了解一下。
线性插值法的通用公式如下：
$$\frac{f(x_2)-f(x_1)}{x_2-x_1}=\frac{f(x)-f(x_1)}{x-x_1}$$
$$f(x)=\frac{x_2-x}{x_2-x_1}f(x_1)+\frac{x-x_1}{x_2-x_1}f(x_2)$$
上面的第二个公式就是权重公式。
所以$$f(0,0.67)=\frac{1-0.67}{1-0}f(0,0)+\frac{0.67-0}{1-0}f(0,1)=1$$
我们再举一个列子：

注意到此时计算出的x为1.3，在1和2之间；y为2.67，在2和3之间。所以我们要通过双线性插值来计算(1.3,2.67)坐标处的像素值。双线性插值算法如下图。

具体就是先对$Q_{11}$、$Q_{21}$进行线性插值，计算出$R_1$，然后对$Q_{12}$、$Q_{22}$进行线性插值，计算出$R_2$，最后对$R_2$、$R_1$进行插值，计算出我们想要的坐标$P$的值。
$$f(1.3,2)=(2-1.3)f(2,2)+(1.3-1)f(1,2)=1$$
$$f(1.3,3)=(2-1.3)f(2,3)+(1.3-1)f(1,3)=1$$
我们对(1.3,2)和(1.3,3)进行线性插值。
$$f(1.3,2.67)=(3-2.67)f(1.3,2)+(2.67-2)f(1.3,3)=1$$
所以根据双线性插值计算出来的(1.3,2.37)处的像素值为1。

总而言之，如果计算出的点在f中是边缘点，则使用单线性插值，如果不是边缘点，则使用双线性插值。
C++代码如下：

int main()
{
    cv::Mat image = cv::imread("LenaRGB.bmp");
    //Using gray image for easy calculations
    cv::Mat grayImage(image.size(), CV_8UC1);
    cv::cvtColor(image, grayImage, CV_BGR2GRAY);
    
    int width = grayImage.cols;
    int height = grayImage.rows;

    // Scale factor
    double k1 = 1.4;
    double k2 = 1.7;
    //Create zero mat to save the outputs
    cv::Mat outImage = cv::Mat::zeros(round(height*k1), round(width*k2), CV_64FC1);

    for (int row = 0; row < outImage.rows; row++)
    {
        for (int col = 0; col < outImage.cols; col++)
        {
        	//得到放大图坐标在原图中对应的坐标
            double srcX = (1 / 1.4)*row;
            double srcY = (1 / 1.7)*col;

            if (srcX > (grayImage.rows - 1))
            {
                srcX = grayImage.rows - 1;
            }
            if (srcY > (grayImage.cols - 1))
            {
                srcY = grayImage.cols - 1;
            }

            // Get x0,x1,y0,y1
            int x0 = floor(srcX);
            int x1 = ceil(srcX);
            int y0 = floor(srcY);
            int y1 = ceil(srcY);

            if ((x0 == x1) && (y0 == y1))
            {
                outImage.at<double>(row, col) = (double)grayImage.at<uchar>(x0, y0);
                continue;
            }
            // 如果是边缘点，则使用单线性插值
            if (x0 == x1)
            {
                double temp = (y1 - srcY)*grayImage.at<uchar>(x0, y0)+
                               (srcY-y0)*grayImage.at<uchar>(x0,y1);
                outImage.at<double>(row, col) = temp;
                continue;
            }
            if (y0 == y1)
            {
                double temp = (x1 - srcX)*grayImage.at<uchar>(x0, y0) +
                    (srcX - x0)*grayImage.at<uchar>(x1, y1);
                outImage.at<double>(row, col) = temp;
                continue;
            }
            // 不是边缘点则使用双线性插值
            double temp1 = (y1 - srcY)*grayImage.at<uchar>(x0, y0) +
                (srcY - y0)*grayImage.at<uchar>(x0, y1);
            double temp2 = (y1 - srcY)*grayImage.at<uchar>(x1, y0) +
                (srcY - y0)*grayImage.at<uchar>(x1, y1);
            double temp = (x1 - srcX)*temp1 + (srcX - x0)*temp2;
            outImage.at<double>(row, col) = temp;
        }
    }
    // Convert CV_64FC1 to CV_8UC1 
    outImage.convertTo(outImage, CV_8UC1);
    return 0;
}