1 算法原理
CBF是*Cross Bilateral Filter(交叉双边滤波)*的缩写,论文《IMAGE FUSION BASED ON PIXEL SIGNIFICANCE USING CROSS BILATERAL FILTER》。
论文中,作者使用交叉双边滤波算法对原始图像
A
A
A,
B
B
B 进行处理得到细节(detail)层,通过细节层得到权重系数,最后经过权重系数对原始图像
A
A
A,
B
B
B进行融合,得到最终的输出图像
F
F
F 。
2 算法实现
2.1 交叉双边滤波算法
双边波算法是局部的,非线性、非迭代的方法,使用高斯滤波器和颜色滤波器来对图像进行处理, 兼具平滑和保持边界的功能。这在我之前的专栏中有提到过,我在这就不在细说双边滤波,接下来说一下文中提到的交叉双边滤波。
交叉双边滤波算法考虑了
A
A
A图像的灰度层的相似性和几何层的闭合性来对
B
B
B 图像进行过滤,或者使用
B
B
B图像来对
A
A
A 图像进行过滤,输入为
A
A
A,
B
B
B 两幅输入图像,假设使用
A
A
A 图像对
B
B
B图像进行过滤,论文中
σ
s
=
1.8
,
σ
r
=
25
\sigma_s=1.8,\sigma_r=25
σs=1.8,σr=25 ,CBFlv滤波窗口大小为
11
∗
11
11*11
11∗11 那么:
B
C
B
F
(
p
)
=
1
W
∑
q
∈
S
(
e
−
∣
p
−
q
∣
2
2
σ
s
2
∗
e
−
∣
A
(
p
)
−
A
(
q
)
∣
2
2
σ
r
2
∗
B
(
q
)
)
B_{CBF}(p)=\frac{1}{W}\sum_{q\in S}(e^{-\frac{|p-q|^2}{2\sigma_s^2}}*e^{-\frac{|A(p)-A(q)|^2}{2\sigma_r^2}}*B(q))
BCBF(p)=W1q∈S∑(e−2σs2∣p−q∣2∗e−2σr2∣A(p)−A(q)∣2∗B(q))
W
=
∑
q
∈
S
e
−
∣
p
−
q
∣
2
2
σ
s
2
∗
e
−
∣
A
(
p
)
−
A
(
q
)
∣
2
2
σ
r
2
W=\sum_{q\in S}e^{-\frac{|p-q|^2}{2\sigma_s^2}}*e^{-\frac{|A(p)-A(q)|^2}{2\sigma_r^2}}
W=q∈S∑e−2σs2∣p−q∣2∗e−2σr2∣A(p)−A(q)∣2
同理,对于
A
C
B
F
(
p
)
A_{CBF}(p)
ACBF(p) 使用的是
B
B
B图像来对
A
A
A 图像进行过滤,由 $B_{CBF}、 A_{CBF} $以及原始图像
A
A
A 、
B
B
B ,在原始图像的基础上减去交叉双边滤波得到的图像就是细节层图像,
A
D
=
A
−
A
C
B
F
,
B
D
=
B
−
B
C
B
F
A_D=A-A_{CBF},B_D=B-B_{CBF}
AD=A−ACBF,BD=B−BCBF
文中作者给出交叉双边滤波有用的解释是:在多聚焦图像中,
A
A
A 图像未聚焦的区域在图像
B
B
B 中是聚焦的区域,对图像
B
B
B 进行交叉双边滤波会使得对图像
B
B
B 的聚焦区域平滑的比图像
B
B
B 的非聚焦区域更强,原因是:图像
A
A
A 中的非聚焦区域由相近的颜色会使得过滤器更像是高斯过滤器。
最上面两幅为原始图像,中间两幅为CBF滤波后的图像,最下面两幅就是由原始图像-CBF图像得到的图像
2.2基于像素的融合策略
接下来的一步是使用细节图像计算融合权重。
选定窗口大小为
w
∗
w
w∗w
w∗w(论文中使用
5
∗
5
5*5
5∗5),来方便计算权重,这个窗口划过的区域看做矩阵 $C_h{i,j}C_h{i,j} $的无偏估计:
c o v a r i a n c e ( X ) = E [ ( X − E [ X ] ) ( X − E ( X ) ) T ] covariance(X)=E[(X-E[X])(X-E(X))^T] covariance(X)=E[(X−E[X])(X−E(X))T]
C h i , j = ∑ k = 1 w ( x k − x ˉ ) ( x k − x ˉ ) T w − 1 C_h^{i,j}=\frac{\sum_{k=1}^w(x_k-\bar{x})(x_k-\bar{x})^T}{w-1} Chi,j=w−1∑k=1w(xk−xˉ)(xk−xˉ)T
得到
C
h
i
,
j
C_h^{i,j}
Chi,j 大小为
w
∗
w
w*w
w∗w,然后计算矩阵
C
h
i
,
j
C_h^{i,j}
Chi,j 的特征值,并且将这些特征值相加就得到了水平方向区域细节强度,并记为
H
d
e
t
a
i
l
S
t
r
e
n
g
t
h
HdetailStrength
HdetailStrength :
H
d
e
t
a
i
l
S
t
r
e
n
g
t
h
(
i
,
j
)
=
∑
k
=
1
w
e
i
g
e
n
k
HdetailStrength(i,j)=\sum_{k=1}^weigen_k
HdetailStrength(i,j)=k=1∑weigenk
同理将 X X X的每一列看做观察量,将 X X X的行看做是变量,计算协方差矩阵 C v i , j C_v^{i,j} Cvi,j,并计算垂直方向区域细节强度 V d e t a i l S t r e n g t h VdetailStrength VdetailStrength :
V d e t a i l S t r e n g t h ( i , j ) = ∑ k = 1 w e i g e n k VdetailStrength(i,j)=\sum_{k=1}^weigen_k VdetailStrength(i,j)=∑k=1weigenk$
最终的权重为,水平和垂直方向的累加和:
w
t
(
i
,
j
)
=
H
d
e
t
a
i
l
S
t
r
e
n
g
t
(
i
,
j
)
+
V
d
e
t
a
i
l
S
t
r
e
n
g
t
(
i
,
j
)
wt(i,j)=HdetailStrengt(i,j)+VdetailStrengt(i,j)
wt(i,j)=HdetailStrengt(i,j)+VdetailStrengt(i,j)
得到 A A A 、 B B B 图像的权重系数后,就可以按照权重系数就行融合:
F ( i , j ) = A ( i , j ) w t a ( i , j ) + B ( i , j ) w t b ( i , j ) w t a ( i , j ) + w t b ( i , j ) F(i,j)=\frac{A(i,j)wt_a(i,j)+B(i,j)wt_b(i,j)}{wt_a(i,j)+wt_b(i,j)} F(i,j)=wta(i,j)+wtb(i,j)A(i,j)wta(i,j)+B(i,j)wtb(i,j)
3 opencv实现代码
import numpy as np
import cv2
import argparse
import math
cov_wsize = 5
sigmas = 1.8
sigmar = 25
ksize = 11
def gaussian_kernel_2d_opencv(kernel_size = 11,sigma = 1.8):
kx = cv2.getGaussianKernel(kernel_size,sigma)
ky = cv2.getGaussianKernel(kernel_size,sigma)
return np.multiply(kx,np.transpose(ky))
def bilateralFilterEx(img_r, img_v):
#edge solved
win_size = ksize//2
img_r_copy = None
img_v_copy = None
img_r_copy = cv2.copyTo(img_r, None)
img_v_copy = cv2.copyTo(img_v, None)
img_r_cbf = np.ones_like(img_r, dtype=np.float32)
img_v_cbf = np.ones_like(img_r, dtype=np.float32)
img_r_copy = np.pad(img_r_copy, (win_size, win_size), 'reflect')
img_v_copy = np.pad(img_v_copy, (win_size, win_size), 'reflect')
gk = gaussian_kernel_2d_opencv()
for i in range(win_size, win_size+img_r.shape[0]):
for j in range(win_size, win_size+img_r.shape[1]):
sumr1 = 0.
sumr2 = 0.
sumv1 = 0.
sumv2 = 0.
img_r_cdis = img_r_copy[i-win_size:i+win_size+1, j-win_size:j+win_size+1] *1.0- img_r_copy[i,j]*1.0
img_v_cdis = img_v_copy[i-win_size:i+win_size+1, j-win_size:j+win_size+1] *1.0- img_v_copy[i,j]*1.0
sumr1 = np.sum(np.exp(-img_v_cdis*img_v_cdis) *gk/ (2*sigmar*sigmar) )
sumv1 = np.sum(np.exp(-img_r_cdis*img_r_cdis) *gk/ (2*sigmar*sigmar) )
sumr2 = np.sum(np.exp(-img_v_cdis*img_v_cdis) *gk*img_r_copy[i-win_size:i+win_size+1, j-win_size:j+win_size+1] *1.0/ (2*sigmar*sigmar) )
sumv2 = np.sum(np.exp(-img_r_cdis*img_r_cdis) *gk*img_v_copy[i-win_size:i+win_size+1, j-win_size:j+win_size+1] *1.0/ (2*sigmar*sigmar) )
img_r_cbf[i-win_size,j-win_size] = sumr2 / sumr1
img_v_cbf[i-win_size,j-win_size] = sumv2 / sumv1
return (img_r*1. - img_r_cbf, img_v*1. - img_v_cbf)
def CBF_WEIGHTS(img_r_d, img_v_d):
win_size = cov_wsize // 2
img_r_weights = np.ones_like(img_r_d, dtype=np.float32)
img_v_weights= np.ones_like(img_v_d, dtype=np.float32)
img_r_d_pad = np.pad(img_r_d, (win_size, win_size), 'reflect')
img_v_d_pad = np.pad(img_v_d, (win_size, win_size), 'reflect')
for i in range(win_size, win_size+img_r_d.shape[0]):
for j in range(win_size, win_size+img_r_d.shape[1]):
npt_r = img_r_d_pad[i-win_size:i+win_size+1, j-win_size:j+win_size+1]
npt_v = img_v_d_pad[i-win_size:i+win_size+1, j-win_size:j+win_size+1]
npt_r_V = npt_r - np.mean(npt_r, axis=0)
npt_r_V = npt_r_V*npt_r_V.transpose()
npt_r_H = npt_r.transpose() - np.mean(npt_r, axis=1)
npt_r_H = npt_r_H*npt_r_H.transpose()
npt_v_V = npt_v - np.mean(npt_v, axis=0)
npt_v_V = npt_v_V*npt_v_V.transpose()
npt_v_H = npt_v.transpose() - np.mean(npt_v, axis=1)
npt_v_H = npt_v_H*npt_v_H.transpose()
img_r_weights[i-win_size,j-win_size] = np.trace(npt_r_H) + np.trace(npt_r_V)
img_v_weights[i-win_size,j-win_size] = np.trace(npt_v_H) + np.trace(npt_v_V)
return img_r_weights, img_v_weights
def CBF_GRAY(img_r, img_v):
img_r_d, img_v_d = bilateralFilterEx(img_r, img_v)
img_r_weights, img_v_weights = CBF_WEIGHTS(img_r_d, img_v_d)
img_fused =(img_r*1. * img_r_weights + img_v*1.*img_v_weights) /(img_r_weights+img_v_weights)
img_fused = cv2.convertScaleAbs(img_fused)
return img_fused
def CBF_RGB(img_r, img_v):
img_r_gray = cv2.cvtColor(img_r, cv2.COLOR_BGR2GRAY)
img_v_gray = cv2.cvtColor(img_v, cv2.COLOR_BGR2GRAY)
return CBF_GRAY(img_r_gray, img_v_gray)
def CBF(_rpath, _vpath):
img_r = cv2.imread(_rpath)
img_v = cv2.imread(_vpath)
if not isinstance(img_r, np.ndarray) :
print('img_r is null')
return
if not isinstance(img_v, np.ndarray) :
print('img_v is null')
return
if img_r.shape[0] != img_v.shape[0] or img_r.shape[1] != img_v.shape[1]:
print('size is not equal')
return
fused_img = None
if len(img_r.shape) < 3 or img_r.shape[2] ==1:
if len(img_v.shape) < 3 or img_v.shape[-1] ==1:
fused_img = CBF_GRAY(img_r, img_v)
else:
img_v_gray = cv2.cvtColor(img_v, cv2.COLOR_BGR2GRAY)
fused_img = CBF_GRAY(img_r, img_v)
else:
if len(img_v.shape) < 3 or img_v.shape[-1] ==1:
img_r_gray = cv2.cvtColor(img_r, cv2.COLOR_BGR2GRAY)
fused_img = CBF_GRAY(img_r_gray, img_v)
else:
fused_img = CBF_RGB(img_r, img_v)
cv2.imshow('fused image', fused_img)
cv2.imwrite("fused_image_2.jpg", fused_img)
cv2.waitKey(0)
if __name__ == '__main__':
parser = argparse.ArgumentParser()
parser.add_argument('-r', type=str, default='ir2.png' ,help='input IR image path', required=False)
parser.add_argument('-v', type=str, default= 'vr2.png',help='input Visible image path', required=False)
args = parser.parse_args()
CBF(args.r, args.v)