Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly
Appearance settings

Commit b9278e4

Browse filesBrowse files
heqin-zhuheqin-zhu
committed
Update image-processing lab: implementation of fft
1 parent fd46d6d commit b9278e4
Copy full SHA for b9278e4

File tree

Expand file treeCollapse file tree

7 files changed

+132
-16
lines changed
Filter options
Expand file treeCollapse file tree

7 files changed

+132
-16
lines changed

‎计算机图像学/labs/dft.py

Copy file name to clipboard
+119Lines changed: 119 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,119 @@
1+
''' mbinary
2+
#########################################################################
3+
# File : fft.py
4+
# Author: mbinary
5+
# Mail: zhuheqin1@gmail.com
6+
# Blog: https://mbinary.xyz
7+
# Github: https://github.com/mbinary
8+
# Created Time: 2019-06-11 12:48
9+
# Description:
10+
#########################################################################
11+
'''
12+
import numpy as np
13+
14+
15+
def _fft_n2(a, invert):
16+
'''O(n^2)'''
17+
N = len(a)
18+
w = np.arange(N)
19+
i = 2j if invert else -2j
20+
m = w.reshape((N, 1)) * w
21+
W = np.exp(m * i * np.pi / N)
22+
return np.concatenate(np.dot(W, a.reshape((N,
23+
1)))) # important, cannot use *
24+
25+
26+
def _fft(a, invert=False):
27+
'''recursion version'''
28+
N = len(a)
29+
if N == 1:
30+
return [a[0]]
31+
elif N & (N - 1) == 0: # O(nlogn), 2^k
32+
even = _fft(a[::2], invert)
33+
odd = _fft(a[1::2], invert)
34+
i = 2j if invert else -2j
35+
factor = np.exp(i * np.pi * np.arange(N // 2) / N)
36+
prod = factor * odd
37+
return np.concatenate([even + prod, even - prod])
38+
else:
39+
return _fft_n2(a, invert)
40+
41+
42+
def _fft2(a, invert=False):
43+
''' iteration version'''
44+
45+
def rev(x):
46+
ret = 0
47+
for i in range(r):
48+
ret <<= 1
49+
if x & 1:
50+
ret += 1
51+
x >>= 1
52+
return ret
53+
54+
N = len(a)
55+
if N & (N - 1) == 0: # O(nlogn), 2^k
56+
r = int(np.log(N))
57+
c = np.array(a,dtype='complex')
58+
i = 2j if invert else -2j
59+
w = np.exp(i * np.pi / N)
60+
for h in range(r - 1, -1, -1):
61+
p = 2**h
62+
z = w**(N / p / 2)
63+
for k in range(N):
64+
if k % p == k % (2 * p):
65+
c[k], c[k + p] = c[k] + c[k + p], c[k] * z**(k % p)
66+
67+
return np.asarray([c[rev(i)] for i in range(N)])
68+
else: # O(n^2)
69+
return _fft_n2(a, invert)
70+
71+
72+
def fft(a):
73+
'''fourier[a]'''
74+
n = len(a)
75+
if n == 0:
76+
raise Exception("[Error]: Invalid length: 0")
77+
return _fft(a)
78+
79+
80+
def ifft(a):
81+
'''invert fourier[a]'''
82+
n = len(a)
83+
if n == 0:
84+
raise Exception("[Error]: Invalid length: 0")
85+
return _fft(a, True) / n
86+
87+
88+
def fft2(arr):
89+
return np.apply_along_axis(fft, 0,
90+
np.apply_along_axis(fft, 1, np.asarray(arr)))
91+
92+
93+
def ifft2(arr):
94+
return np.apply_along_axis(ifft, 0,
95+
np.apply_along_axis(ifft, 1, np.asarray(arr)))
96+
97+
98+
def test(n=128):
99+
print('\nsequence length:', n)
100+
print('fft')
101+
li = np.random.random(n)
102+
print(np.allclose(fft(li), np.fft.fft(li)))
103+
104+
print('ifft')
105+
li = np.random.random(n)
106+
print(np.allclose(ifft(li), np.fft.ifft(li)))
107+
108+
print('fft2')
109+
li = np.random.random(n * n).reshape((n, n))
110+
print(np.allclose(fft2(li), np.fft.fft2(li)))
111+
112+
print('ifft2')
113+
li = np.random.random(n * n).reshape((n, n))
114+
print(np.allclose(ifft2(li), np.fft.ifft2(li)))
115+
116+
117+
if __name__ == '__main__':
118+
for i in range(1, 4):
119+
test(i * 16)

‎计算机图像学/labs/lab1.py

Copy file name to clipboardExpand all lines: 计算机图像学/labs/lab1.py
+1-1Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -103,7 +103,7 @@ def show(img, s='opencv'):
103103

104104
cmap = mpl.cm.gray # mpl.cm.gray_r 'gray'
105105

106-
plt.figure(figsize=(10, 10))
106+
plt.figure(figsize=(8, 8))
107107

108108
plt.subplot(321), plt.imshow(img,cmap=cmap), plt.title('origin'),plt.xticks([]), plt.yticks([])
109109
plt.subplot(322), plt.imshow(img2,cmap=cmap), plt.title(f'tran k={k},b={b}'),plt.xticks([]), plt.yticks([])

‎计算机图像学/labs/lab2.py

Copy file name to clipboardExpand all lines: 计算机图像学/labs/lab2.py
+1-1Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -142,7 +142,7 @@ def find_median(arr):
142142
img3 = median_filter(noised_img)
143143

144144
cmap = mpl.cm.gray
145-
plt.figure(figsize=(10, 10))
145+
plt.figure(figsize=(8, 8))
146146

147147
plt.subplot(221), plt.xticks([]), plt.yticks([])
148148
plt.imshow(img, cmap=cmap)

‎计算机图像学/labs/lab3.py

Copy file name to clipboardExpand all lines: 计算机图像学/labs/lab3.py
+1-1Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -46,7 +46,7 @@ def prewitt(img):
4646
img2 = roberts(img)
4747
img3 = prewitt(img)
4848
cmap = mpl.cm.gray
49-
plt.figure(figsize=(10, 10))
49+
plt.figure(figsize=(8, 8))
5050

5151
plt.subplot(221),plt.xticks([]), plt.yticks([])
5252
plt.imshow(img, cmap=cmap)

‎计算机图像学/labs/lab4.py

Copy file name to clipboardExpand all lines: 计算机图像学/labs/lab4.py
+10-13Lines changed: 10 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -4,36 +4,33 @@
44
import numpy as np
55
import matplotlib as mpl
66
from matplotlib import pyplot as plt
7+
import dft
78

89

9-
def fft(img):
10-
f = np.fft.fft2(img)
11-
# fshift = np.fft.fftshift(f)
12-
# magnitude_spectrum
13-
mag = np.abs(f)
14-
imag = np.fft.ifft2(mag)
15-
phase = np.angle(f)
16-
iphase = np.fft.ifft2(phase)
17-
return mag, phase, imag, iphase
1810

1911

2012
if __name__ == '__main__':
2113
path = sys.argv[1]
2214
img = cv2.imread(path, cv2.IMREAD_GRAYSCALE)
2315

24-
mag, phase, imag, iphase = fft(img)
16+
f = dft.fft2(img)
17+
# magnitude_spectrum
18+
mag = np.abs(f)
19+
invert_magnitude = dft.ifft2(mag)
2520
mag = np.log(mag + 1)
26-
imag = np.real(imag)
21+
phase = np.angle(f)
22+
iphase = dft.ifft2(phase)
23+
invert_magnitude = np.real(invert_magnitude)
2724
iphase = np.real(iphase)
2825

2926
cmap = mpl.cm.gray
30-
plt.figure(figsize=(10, 10))
27+
plt.figure(figsize=(6,6))
3128

3229
plt.subplot(221), plt.imshow(img, cmap=cmap)
3330
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
3431
plt.subplot(222), plt.imshow(mag, cmap=cmap)
3532
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
36-
plt.subplot(223), plt.imshow(imag, cmap=cmap)
33+
plt.subplot(223), plt.imshow(invert_magnitude, cmap=cmap)
3734
plt.title('idft [ magnitude] '), plt.xticks([]), plt.yticks([])
3835
plt.subplot(224), plt.imshow(iphase, cmap=cmap)
3936
plt.title('idft [ phase ]'), plt.xticks([]), plt.yticks([])

‎计算机图像学/labs/result/lab4-rect1.png

Copy file name to clipboard
-17.9 KB
Loading

‎计算机图像学/labs/result/lab4-rect2.png

Copy file name to clipboard
-25 KB
Loading

0 commit comments

Comments
0 (0)
Morty Proxy This is a proxified and sanitized view of the page, visit original site.