I thought I'd better try this sparse Gabor filtering approximation algorithm before I got too proud of it. It does seem to more or less work in this case, but I haven't calculated its error or anything.

In [1]:
%pylab inline

n=5
m=3
p=30
Populating the interactive namespace from numpy and matplotlib
In [25]:
impulse = zeros((700, 700), dtype=dtype('float16'))
impulse[300, 100] = 1.0
ir1 = impulse.copy()
for i in range(4*n, len(ir1)):
    ir1[i, 3*n:] -= ir1[i-4*n, :-3*n]
contour(ir1); colorbar()
matshow(ir1[290:350, 90:150], origin='lower'); colorbar()
ir1.min(), ir1.max()
Out[25]:
(-1.0, 1.0)
In [3]:
fr1 = ir1.copy()
for i in range(8*n*m, len(fr1)):
    fr1[i, 6*n*m:] -= ir1[i-8*n*m, :-6*n*m]
contour(fr1)
(abs(fr1) > 0).sum() # should be 2*m
Out[3]:
6
In [4]:
ir2 = fr1.copy()
for i in range(4*n, len(ir2)):
    ir2[i, 3*n:] -= ir2[i-4*n, :-3*n]
fr2 = ir2.copy()
for i in range(8*n*m, len(fr2)):
    fr2[i, 6*n*m:] -= ir2[i-8*n*m, :-6*n*m]
contour(fr2); colorbar()
(abs(fr2) > 0).sum() # should be 4*m-1 = 11?
Out[4]:
11
In [5]:
ir3 = fr2.copy()
for i in range(4*n, len(ir3)):
    ir3[i, 3*n:] -= ir3[i-4*n, :-3*n]
fr3 = ir3.copy()
for i in range(8*n*m, len(fr3)):
    fr3[i, 6*n*m:] -= ir3[i-8*n*m, :-6*n*m]
contour(fr3)
(abs(fr3) > 0).sum() # should be 6*m-2 = 16?
Out[5]:
16
In [6]:
plot(fr3.sum(axis=1))
Out[6]:
[<matplotlib.lines.Line2D at 0xb081690c>]
In [7]:
%%time
hs = fr3.copy()
for i in range(int(4*n/3), len(hs)):
    hs[i, 3*n/3:] += fr3[i-int(4*n/3), :-3*n/3]
subplot(121); contour(hs)
subplot(122); plot(hs.sum(axis=1))
CPU times: user 1.32 s, sys: 68 ms, total: 1.38 s
Wall time: 1.54 s
In [8]:
contour(hs[300:700, :400]); colorbar()
Out[8]:
<matplotlib.colorbar.Colorbar instance at 0xb0803b2c>
In [9]:
ai1 = hs.copy()
for i in range(4, len(ai1)):
    ai1[i, 3:] += ai1[i-4, :-3]
ff1 = ai1.copy()
ff1[8:, 6:] -= ai1[:-8, :-6]
contour(ff1)
Out[9]:
<matplotlib.contour.QuadContourSet instance at 0xb0fa4cec>
In [10]:
contour(ff1[200:400, :200])
colorbar()
Out[10]:
<matplotlib.colorbar.Colorbar instance at 0xb09c81ac>
In [11]:
ai2 = ff1.copy()
for i in range(4, len(ai2)):
    ai2[i, 3:] += ai2[i-4, :-3]
ff2 = ai2.copy()
ff2[8:, 6:] -= ai2[:-8, :-6]
contour(ff2); colorbar(); ff2.min(), ff2.max(), ff2.sum()
Out[11]:
(-54.0, 54.0, 0.0)
In [12]:
plot([ff2r.sum() for ff2r in ff2[200:700]], '.')
Out[12]:
[<matplotlib.lines.Line2D at 0xaf961aac>]
In [13]:
bf1i = ff2.copy()
for j in range(1, len(bf1i[0])):
    bf1i[:, j] += bf1i[:, j-1]
bf1c = bf1i.copy()
bf1c[:, 8:] -= bf1i[:, :-8]
contour(bf1c)
Out[13]:
<matplotlib.contour.QuadContourSet instance at 0xb09f016c>
In [14]:
bf2i = bf1c.copy()
for j in range(1, len(bf2i[0])):
    bf2i[:, j] += bf2i[:, j-1]
bf2c = bf2i.copy()
bf2c[:, 8:] -= bf2i[:, :-8]
contour(bf2c)
Out[14]:
<matplotlib.contour.QuadContourSet instance at 0xb08dbd4c>
In [15]:
# now vertical
bf3i = bf2c.copy()
for i in range(1, len(bf3i)):
    bf3i[i] += bf3i[i-1]
bf3c = bf3i.copy()
bf3c[6:] -= bf3i[:-6]
contour(bf3c)
Out[15]:
<matplotlib.contour.QuadContourSet instance at 0xaf7748cc>
In [16]:
bf4i = bf3c.copy()
for i in range(1, len(bf4i)):
    bf4i[i] += bf4i[i-1]
bf4c = bf4i.copy()
bf4c[6:] -= bf4i[:-6]
contour(bf4c)
Out[16]:
<matplotlib.contour.QuadContourSet instance at 0xadb877ec>
In [17]:
gca().set_aspect('equal')
contour(bf4c[300:500, 100:300])
colorbar()
Out[17]:
<matplotlib.colorbar.Colorbar instance at 0xadb115cc>
In [18]:
plot([bf4cr.sum() for bf4cr in bf4c[300:600]], '.')
Out[18]:
[<matplotlib.lines.Line2D at 0xacd4572c>]
In [19]:
tvi1 = bf4c.copy()
for j in range(4, len(tvi1[0])):
    tvi1[:-3, j] += tvi1[3:, j-4]
tvc1 = tvi1.copy()
tvc1[:-3*p, 4*p:] -= tvi1[3*p:, :-4*p]
gca().set_aspect('equal')
contour(tvc1); colorbar(); tvc1.max(), tvc1.min()
Out[19]:
(6216.0, -6216.0)
In [20]:
tvi2 = tvc1 / 256
for j in range(4, len(tvi2[0])):
    tvi2[:-3, j] += tvi2[3:, j-4]
tvc2 = tvi2.copy()
tvc2[:-3*p, 4*p:] -= tvi2[3*p:, :-4*p]
gca().set_aspect('equal')
contour(tvc2)
Out[20]:
<matplotlib.contour.QuadContourSet instance at 0xabe171cc>
In [21]:
tvi3 = tvc2 / 256
for j in range(4, len(tvi3[0])):
    tvi3[:-3, j] += tvi3[3:, j-4]
tvc3 = tvi3.copy()
tvc3[:-3*p, 4*p:] -= tvi3[3*p:, :-4*p]
gca().set_aspect('equal')
contour(tvc3)
Out[21]:
<matplotlib.contour.QuadContourSet instance at 0xabf7756c>
In [22]:
gca().set_aspect('equal')
contour(tvc3[200:500, 300:600])
colorbar()
tvc3.min(), tvc3.max()
Out[22]:
(-64.003784, 64.003784)
In [24]:
imshow(tvc3[200:500, 300:600], origin='lower'); colorbar()
Out[24]:
<matplotlib.colorbar.Colorbar instance at 0xab2078ec>
In [23]: