• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1#!/usr/bin/env python
2
3'''
4Wiener deconvolution.
5
6Sample shows how DFT can be used to perform Weiner deconvolution [1]
7of an image with user-defined point spread function (PSF)
8
9Usage:
10  deconvolution.py  [--circle]
11      [--angle <degrees>]
12      [--d <diameter>]
13      [--snr <signal/noise ratio in db>]
14      [<input image>]
15
16  Use sliders to adjust PSF paramitiers.
17  Keys:
18    SPACE - switch btw linear/cirular PSF
19    ESC   - exit
20
21Examples:
22  deconvolution.py --angle 135 --d 22  ../data/licenseplate_motion.jpg
23    (image source: http://www.topazlabs.com/infocus/_images/licenseplate_compare.jpg)
24
25  deconvolution.py --angle 86 --d 31  ../data/text_motion.jpg
26  deconvolution.py --circle --d 19  ../data/text_defocus.jpg
27    (image source: compact digital photo camera, no artificial distortion)
28
29
30[1] http://en.wikipedia.org/wiki/Wiener_deconvolution
31'''
32
33import numpy as np
34import cv2
35
36# local module
37from common import nothing
38
39
40def blur_edge(img, d=31):
41    h, w  = img.shape[:2]
42    img_pad = cv2.copyMakeBorder(img, d, d, d, d, cv2.BORDER_WRAP)
43    img_blur = cv2.GaussianBlur(img_pad, (2*d+1, 2*d+1), -1)[d:-d,d:-d]
44    y, x = np.indices((h, w))
45    dist = np.dstack([x, w-x-1, y, h-y-1]).min(-1)
46    w = np.minimum(np.float32(dist)/d, 1.0)
47    return img*w + img_blur*(1-w)
48
49def motion_kernel(angle, d, sz=65):
50    kern = np.ones((1, d), np.float32)
51    c, s = np.cos(angle), np.sin(angle)
52    A = np.float32([[c, -s, 0], [s, c, 0]])
53    sz2 = sz // 2
54    A[:,2] = (sz2, sz2) - np.dot(A[:,:2], ((d-1)*0.5, 0))
55    kern = cv2.warpAffine(kern, A, (sz, sz), flags=cv2.INTER_CUBIC)
56    return kern
57
58def defocus_kernel(d, sz=65):
59    kern = np.zeros((sz, sz), np.uint8)
60    cv2.circle(kern, (sz, sz), d, 255, -1, cv2.LINE_AA, shift=1)
61    kern = np.float32(kern) / 255.0
62    return kern
63
64
65if __name__ == '__main__':
66    print __doc__
67    import sys, getopt
68    opts, args = getopt.getopt(sys.argv[1:], '', ['circle', 'angle=', 'd=', 'snr='])
69    opts = dict(opts)
70    try:
71        fn = args[0]
72    except:
73        fn = '../data/licenseplate_motion.jpg'
74
75    win = 'deconvolution'
76
77    img = cv2.imread(fn, 0)
78    if img is None:
79        print 'Failed to load fn1:', fn1
80        sys.exit(1)
81
82    img = np.float32(img)/255.0
83    cv2.imshow('input', img)
84
85    img = blur_edge(img)
86    IMG = cv2.dft(img, flags=cv2.DFT_COMPLEX_OUTPUT)
87
88    defocus = '--circle' in opts
89
90    def update(_):
91        ang = np.deg2rad( cv2.getTrackbarPos('angle', win) )
92        d = cv2.getTrackbarPos('d', win)
93        noise = 10**(-0.1*cv2.getTrackbarPos('SNR (db)', win))
94
95        if defocus:
96            psf = defocus_kernel(d)
97        else:
98            psf = motion_kernel(ang, d)
99        cv2.imshow('psf', psf)
100
101        psf /= psf.sum()
102        psf_pad = np.zeros_like(img)
103        kh, kw = psf.shape
104        psf_pad[:kh, :kw] = psf
105        PSF = cv2.dft(psf_pad, flags=cv2.DFT_COMPLEX_OUTPUT, nonzeroRows = kh)
106        PSF2 = (PSF**2).sum(-1)
107        iPSF = PSF / (PSF2 + noise)[...,np.newaxis]
108        RES = cv2.mulSpectrums(IMG, iPSF, 0)
109        res = cv2.idft(RES, flags=cv2.DFT_SCALE | cv2.DFT_REAL_OUTPUT )
110        res = np.roll(res, -kh//2, 0)
111        res = np.roll(res, -kw//2, 1)
112        cv2.imshow(win, res)
113
114    cv2.namedWindow(win)
115    cv2.namedWindow('psf', 0)
116    cv2.createTrackbar('angle', win, int(opts.get('--angle', 135)), 180, update)
117    cv2.createTrackbar('d', win, int(opts.get('--d', 22)), 50, update)
118    cv2.createTrackbar('SNR (db)', win, int(opts.get('--snr', 25)), 50, update)
119    update(None)
120
121    while True:
122        ch = cv2.waitKey() & 0xFF
123        if ch == 27:
124            break
125        if ch == ord(' '):
126            defocus = not defocus
127            update(None)
128