isp算法学习-插值
800 Words|Read in about 4 Min|本文总阅读量次
本文是isp算法学习的第四篇章，主要是讲述插值原理和算法。
1插值方法
1.1邻域插值
这个时候落点P不在区域四个顶点上，那么求其P距离区域顶点的值 $$ min = {r_1, r_2, r_3,r_4} $$ 把整个区域分成几个小区域，然后判断点落在那个区域，就意味着P点离对应的点最近，然后就令P点等于最近的点。
上图中P离区域顶点A最近，所以选择A点作为P的值。
通常在代码中不会使用求距离公式，比较复杂 $$ d = \sqrt{(x-x_0)^2 + (y-y_0)^2} $$ 首先对整个区域归一化处理，然后在matlab会使用建坐标系来简化运算量和步骤。
以区域中心建立xy坐标系，P点通过在四个象限落点来判断P点最近的点。
1.2线性插值
线性插值，就是利用两个点，进行插值，得到一个近似的点，让P等于这个点。
以上图为例，可以选择用A1和A2点来插值，也可以使用A1和A3点来插值，也可以用A1和A4点来插值。
相对于邻域插值来演，用两个点插值会更加精确。
1.3双线性插值
双线性就是，先进行两次线性插值，然后再利先前两次插值的结果再进行一次线性插值，最终得到P
以上图为例，第一次插值，选择用A1和A2点来插值求得B1，第二次插值，选择用A3和A4点来插值求得B2。第三次插值，选择用B1和B2点来插值求得P
相对于单次线性插值来演，双线性插值用四个点插值会更加精确。
进一步的，可以列出公式，归一化之后以左下角为原点建立xyz坐标
$$
\begin{cases}
\frac{Z_b-Z_a}{Y_1-Y_0} = \frac{Z_n-Z_a}{Y-Y_0}&(1)\
\frac{Z_c-Z_d}{Y_1-Y_0} = \frac{Z_m-Z_a}{Y-Y_0}&(2)\
\frac{Z_m-Z_n}{X_1-X_0} = \frac{Z-Z_0}{X-X_0}&(3)
\end{cases}
$$
得到最终z的值
$$
Z = (Z_d-Z_a)* X + (Z_b -Z_a)*Y+(Z_c+Z_a-Z_b-Z_d)*X*Y + Z_a
$$
1.4双三次插值
因为连续信号在采样的过程中会变成离散信号，离散信号变成连续信号的公式，如下所示
$$
x(t) = \sum_{i=-\infty}^{+\infty}{x(nT_s)sinc(
\frac{\pi}{T_s}(t-nT_s)
)}
$$
转变成连续信号后，对连续信号处理，sin函数比较复杂，取上述公式的[-2,2]区间对应的值，如右图所示，那么简化分成三段
$$
S(x) =
\begin{cases}
1-2|x|^2+|x|^3 &|x|<1\
1-2|x|^2+|x|^3 &1\leq|x|\leq2\
0 &|x|>2
\end{cases}
$$
2算法
2.1邻域插值
1clc,clear,close all;
2% 读取图片
3orgImage = imread('./images/lena.bmp');
4figure;imshow(orgImage);title('org image');
5
6% 获取长宽
7[width, height] = size(orgImage);
8m = width / 2;
9n = height / 2;
10smallImage = zeros(m,n);
11% 降采样，将原图缩减为原来的1/2
12for i=1:m
13 for j=1:n
14 smallImage(i,j) = orgImage(2*i,2*j);
15 end
16end
17figure;imshow(uint8(smallImage));title('small image');
18
19% 插值时需要特殊处理四周最外圈的行和列，本算法中将其向外扩展一圈，用最外圈的值填充
20headRowMat = smallImage(1,:);%取f的第1行
21tailRowMat = smallImage(m,:);%取f的第m行
22% 行扩展后，列扩展时需要注意四个角需要单独扩展进去，不然就成了十字架形的
23headColumnMat = [smallImage(1,1), smallImage(:,1)', smallImage(m,1)];
24tailColumnMat = [smallImage(1,n), smallImage(:,n)', smallImage(m,n)];
25expandImage = [headRowMat; smallImage; tailRowMat];
26expandImage = [headColumnMat; expandImage'; tailColumnMat];
27expandImage = uint8(expandImage');
28figure;imshow(expandImage);title('expand image');
29
30% 按比例放大
31[smallWidth, smallHeight] = size(smallImage);
32% 设置放大系数
33magnification = 2;
34newWidth = magnification * smallWidth;
35newHeight = magnification * smallHeight;
36% 创建一个新的矩阵，用于承接变换后的图像
37newImage = zeros(newWidth, newHeight);
38
39% 循环计算出新图像的像素值
40for i = 1 : newWidth
41 for j = 1: newHeight
42 detaX = rem(i, magnification) / magnification;
43 floorX = floor(i / magnification) + 1;
44 detaY = rem(j, magnification) / magnification;
45 floorY = floor(j / magnification) + 1;
46 if detaX < 0.5 && detaY < 0.5
47 newImage(i, j) = expandImage(floorX, floorY);
48 elseif detaX < 0.5 && detaY >= 0.5
49 newImage(i, j) = expandImage(floorX, floorY + 1);
50 elseif detaX >= 0.5 && detaY < 0.5
51 newImage(i, j) = expandImage(floorX + 1, floorY);
52 else
53 newImage(i, j) = expandImage(floorX + 1, floorY + 1);
54 end
55 end
56end
57figure;imshow(uint8(newImage));title('NeighborhoodInterpolation');
58figure;
59subplot(121);imshow(uint8(smallImage));title('small img');
60subplot(122);imshow(uint8(newImage));title('NeighborhoodInterpolation img');
2.2双线性插值
1clc,clear,close all;
2% 读取图片
3orgImage = imread('./images/lena.bmp');
4figure;imshow(orgImage);title('org image');
5
6% 获取长宽
7[width, height] = size(orgImage);
8m = width / 2;
9n = height / 2;
10smallImage = zeros(m,n);
11% 降采样，将原图缩减为原来的1/2
12for i=1:m
13 for j=1:n
14 smallImage(i,j) = orgImage(2*i,2*j);
15 end
16end
17figure;imshow(uint8(smallImage));title('small image');
18
19% 插值时需要特殊处理四周最外圈的行和列，本算法中将其向外扩展一圈，用最外圈的值填充
20headRowMat = smallImage(1,:);%取f的第1行
21tailRowMat = smallImage(m,:);%取f的第m行
22% 行扩展后，列扩展时需要注意四个角需要单独扩展进去，不然就成了十字架形的
23headColumnMat = [smallImage(1,1), smallImage(:,1)', smallImage(m,1)];
24tailColumnMat = [smallImage(1,n), smallImage(:,n)', smallImage(m,n)];
25expandImage = [headRowMat; smallImage; tailRowMat];
26expandImage = [headColumnMat; expandImage'; tailColumnMat];
27expandImage = uint8(expandImage');
28figure;imshow(expandImage);title('expand image');
29
30% 按比例放大
31[smallWidth, smallHeight] = size(smallImage);
32% 设置放大系数
33magnification = 2;
34newWidth = magnification * smallWidth;
35newHeight = magnification * smallHeight;
36% 创建一个新的矩阵，用于承接变换后的图像
37newImage = zeros(newWidth, newHeight);
38
39% f(x,y) = [f(1,0)-f(0,0)]*x+[f(0,1)-f(0,0)]*y+[f(1,1)+f(0,0)-f(1,0)-f(0,1)]*xy+f(0,0)
40for i = 1 : newWidth
41 for j = 1: newHeight
42 detaX = rem(i, magnification) / magnification;
43 floorX = floor(i / magnification) + 1;
44 detaY = rem(j, magnification) / magnification;
45 floorY = floor(j / magnification) + 1;
46 newImage(i, j) = (expandImage(floorX + 1,floorY) - expandImage(floorX,floorY)) * detaX + ...
47 (expandImage(floorX, floorY + 1) - expandImage(floorX, floorY)) * detaY + ...
48 (expandImage(floorX+1, floorY+1) + expandImage(floorX, floorY) - ...
49 expandImage(floorX+1, floorY) - expandImage(floorX, floorY+1)) * detaX * detaY + ...
50 expandImage(floorX, floorY);
51 end
52end
53figure;
54% subplot(121);imshow(uint8(smallImage));title('small image');
55imshow(uint8(newImage));title('BilinearInterpolation');
2.3双三次内插
1clc,clear,close all;
2orgImage = imread('./images/lena.bmp');
3[width, height] = size(orgImage);%将图像隔行隔列抽取元素，得到缩小的图像f
4figure; imshow(orgImage); title('org image');%显示原图像
5
6m = width/2;
7n = height/2;
8smallImage = zeros(m, n);
9for i = 1: m
10 for j = 1: n
11 smallImage(i, j) = orgImage(2*i, 2*j);
12 end
13end
14figure;imshow(uint8(smallImage));title('small image');%显示缩小的图像
15
16
17magnification = 2;%设置放大倍数
18a = smallImage(1,:);%取f的第1行
19c = smallImage(m,:);%取f的第m行
20%将待插值图像矩阵前后各扩展两行两列,共扩展四行四列到f1
21b = [smallImage(1,1), smallImage(1,1), smallImage(:,1)', smallImage(m,1), smallImage(m,1)];
22d = [smallImage(1,n), smallImage(1,n), smallImage(:,n)', smallImage(m,n), smallImage(m,n)];
23a1 = [a; a; smallImage; c; c];
24b1 = [b; b; a1'; d; d];
25expandImage = double(b1');
26
27newImage = zeros(magnification*m,magnification*n);
28for i = 1:magnification * m%利用双三次插值公式对新图象所有像素赋值
29 u = rem(i, magnification)/magnification;
30 i1 = floor(i/magnification) + 2;%floor()向左取整，floor(1.3)=floor(1.7)=1
31 A = [sw(1+u) sw(u) sw(1-u) sw(2-u)];
32 for j = 1:magnification*n
33 v = rem(j, magnification)/magnification; j1=floor(j/magnification)+2;
34 C = [sw(1+v); sw(v); sw(1-v); sw(2-v)];
35 B = [expandImage(i1-1,j1-1) expandImage(i1-1,j1) expandImage(i1-1,j1+1) expandImage(i1-1,j1+2);
36 expandImage(i1,j1-1) expandImage(i1,j1) expandImage(i1,j1+1) expandImage(i1,j1+2);
37 expandImage(i1+1,j1-1) expandImage(i1+1,j1) expandImage(i1+1,j1+1) expandImage(i1+1,j1+2);
38 expandImage(i1+2,j1-1) expandImage(i1+2,j1) expandImage(i1+2,j1+1) expandImage(i1+2,j1+2)];
39 newImage(i,j) = (A*B*C);
40 end
41end
42%显示插值后的图像
43figure,
44% subplot(121);imshow(uint8(smallImage));title('small image');%显示缩小的图像
45imshow(uint8(newImage));title('BicubicInterpolation');
其中函数sw是双三次插值算法sin函数的拟合函数
1function A = sw(w1)
2w = abs(w1);
3if w < 1 && w >= 0
4 A = 1 - 2 * w^2 + w^3;
5elseif w >= 1 && w < 2
6 A = 4 - 8 * w + 5 * w^2 - w^3;
7else
8 A = 0;
9end