-
Notifications
You must be signed in to change notification settings - Fork 3
/
guidedFilter.py
96 lines (78 loc) · 3.33 KB
/
guidedFilter.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
from itertools import combinations_with_replacement
from collections import defaultdict
import numpy as np
from numpy.linalg import inv
B, G, R = 0, 1, 2 # index for convenience
def boxfilter(I, r):
"""Fast box filter implementation.
Parameters
----------
I: a single channel/gray image data normalized to [0.0, 1.0]
r: window radius
Return
-----------
The filtered image data.
"""
M, N = I.shape
dest = np.zeros((M, N))
#print(I)
# cumulative sum over Y axis (tate-houkou no wa)
sumY = np.cumsum(I, axis=0)
#print('sumY:{}'.format(sumY))
# difference over Y axis
dest[:r + 1] = sumY[r:2*r + 1] # top r+1 lines
dest[r + 1:M - r] = sumY[2*r + 1:] - sumY[:M - 2*r - 1]
#print(sumY[2*r + 1:]) # from 2*r+1 to end lines
#print(sumY[:M - 2*r - 1]) # same lines of above, from start
#tile replicate sumY[-1] and line them up to match the shape of (r, 1)
dest[-r:] = np.tile(sumY[-1], (r, 1)) - sumY[M - 2*r - 1:M - r - 1] # bottom r lines
# cumulative sum over X axis
sumX = np.cumsum(dest, axis=1)
#print('sumX:{}'.format(sumX))
# difference over X axis
dest[:, :r + 1] = sumX[:, r:2*r + 1] # left r+1 columns
dest[:, r + 1:N - r] = sumX[:, 2*r + 1:] - sumX[:, :N - 2*r - 1]
dest[:, -r:] = np.tile(sumX[:, -1][:, None], (1, r)) - sumX[:, N - 2*r - 1:N - r - 1] # right r columns
#print(dest)
return dest
def guided_filter(I, p, r=15, eps=1e-3):
"""Refine a filter under the guidance of another (RGB) image.
Parameters
-----------
I: an M * N * 3 RGB image for guidance.
p: the M * N filter to be guided. transmission is used for this case.
r: the radius of the guidance
eps: epsilon for the guided filter
Return
-----------
The guided filter.
"""
M, N = p.shape
base = boxfilter(np.ones((M, N)), r) # this is needed for regularization
# each channel of I filtered with the mean filter. this is myu.
means = [boxfilter(I[:, :, i], r) / base for i in range(3)]
# p filtered with the mean filter
mean_p = boxfilter(p, r) / base
# filter I with p then filter it with the mean filter
means_IP = [boxfilter(I[:, :, i]*p, r) / base for i in range(3)]
# covariance of (I, p) in each local patch
covIP = [means_IP[i] - means[i]*mean_p for i in range(3)]
# variance of I in each local patch: the matrix Sigma in ECCV10 eq.14
var = defaultdict(dict)
for i, j in combinations_with_replacement(range(3), 2):
var[i][j] = boxfilter(I[:, :, i]*I[:, :, j], r) / base - means[i]*means[j]
a = np.zeros((M, N, 3))
for y, x in np.ndindex(M, N):
# rr, rg, rb
# Sigma = rg, gg, gb
# rb, gb, bb
Sigma = np.array([[var[B][B][y, x], var[B][G][y, x], var[B][R][y, x]],
[var[B][G][y, x], var[G][G][y, x], var[G][R][y, x]],
[var[B][R][y, x], var[G][R][y, x], var[R][R][y, x]]])
cov = np.array([c[y, x] for c in covIP])
a[y, x] = np.dot(cov, inv(Sigma + eps*np.eye(3))) # eq 14
# ECCV10 eq.15
b = mean_p - a[:, :, R]*means[R] - a[:, :, G]*means[G] - a[:, :, B]*means[B]
# ECCV10 eq.16
q = (boxfilter(a[:, :, R], r)*I[:, :, R] + boxfilter(a[:, :, G], r)*I[:, :, G] + boxfilter(a[:, :, B], r)*I[:, :, B] + boxfilter(b, r)) / base
return q