12 # size of block that will be matched
13 # (correlate a block this size against a block with added slippage all around)
17 # pitch between the block centres
21 # allowable +/- slippage between pairs
32 def calc_bandpass(image):
33 _, _, cs = image.shape
36 scipy.ndimage.gaussian_filter(image[:, :, i], CUTOFF1, mode = 'mirror') -
37 scipy.ndimage.gaussian_filter(image[:, :, i], CUTOFF0, mode = 'mirror')
43 def calc_match(bandpass0, bandpass1, xc, yc):
46 x1 = xc - XM // 2 - XS
47 y1 = yc - YM // 2 - YS
60 # note: swapping block1, block0 flips the output (subtracts x and y
61 # from BLOCK_SIZE) and we need this to find matching part of image1
65 scipy.signal.correlate(
76 #temp = corr - numpy.mean(corr)
77 #temp /= 10. * numpy.sqrt(numpy.mean(numpy.square(temp)))
78 #gamma.write_image(f'corr_{i:d}_{j:d}.jpg', temp + .5)
80 # find slippage from correlation
81 yo, xo = numpy.unravel_index(numpy.argmax(corr), corr.shape)
82 max_corr = corr[yo, xo]
86 xo > XS * 2 - CUTOFF1 or
92 # estimate position within block of feature being matched
93 block1 = block1[yo:yo + YM, xo:xo + XM, :]
94 match = numpy.sum(block0 * block1, 2)
95 #print('xxx', numpy.sum(match), max_corr)
98 numpy.sum(match, 0) @ numpy.arange(XM, dtype = numpy.double)
101 numpy.sum(match, 1) @ numpy.arange(YM, dtype = numpy.double)
104 # x = int(math.floor(xf))
105 # y = int(math.floor(yf))
107 # diag0 = block0 + .5
109 # max(y - 21, 0):max(y + 22, 0),
110 # max(x - 1, 0):max(x + 2, 0),
114 # max(y - 1, 0):max(y + 2, 0),
115 # max(x - 21, 0):max(x + 22, 0),
118 # gamma.write_image(f'diag_{xc:d}_{yc:d}_0.jpg', diag0)
120 # diag1 = block1 + .5
122 # max(y - 21, 0):max(y + 22, 0),
123 # max(x - 1, 0):max(x + 2, 0),
127 # max(y - 1, 0):max(y + 2, 0),
128 # max(x - 21, 0):max(x + 22, 0),
131 # gamma.write_image(f'diag_{xc:d}_{yc:d}_1.jpg', diag1)
133 # return offset and feature relative to block centre
134 return xo - XS, yo - YS, xf - XM // 2, yf - YM // 2
139 if len(sys.argv) >= 2 and sys.argv[1] == '--diag':
143 in_jpg0 = 'tank_battle/down_2364.jpg'
144 in_jpg1 = 'tank_battle/down_2365.jpg'
145 out_jpg0 = 'out0_2364.jpg'
146 out_jpg1 = 'out0_2365.jpg'
148 print(f'read {in_jpg0:s}')
149 image0 = gamma.read_image(in_jpg0)
153 bandpass0 = calc_bandpass(image0)
155 gamma.write_image('bandpass0.jpg', bandpass0 + .5)
157 print(f'read {in_jpg1:s}')
158 image1 = gamma.read_image(in_jpg1)
159 assert image1.shape == shape
162 bandpass1 = calc_bandpass(image1)
164 gamma.write_image('bandpass1.jpg', bandpass1 + .5)
167 xb = (xs // 2 - XM - 2 * XS) // XP
168 yb = (ys // 2 - YM - 2 * YS) // YP
169 print('xb', xb, 'yb', yb)
171 print('find corner candidates')
174 corner_candidates = []
177 print('i', i, 'j', j)
179 # correlate blocks in (i, j)-corner
183 yc = YS + YM // 2 + k * YP
187 xc = XS + XM // 2 + l * XP
190 match = calc_match(bandpass0, bandpass1, xc, yc)
191 if match is not None:
192 offsets.append(match [:2])
193 xo, yo, xf, yf = match
198 p_all.append(numpy.array([xf0, yf0], numpy.double))
199 q_all.append(numpy.array([xf1, yf1], numpy.double))
200 blocks.append((xf0, yf0, xf1, yf1))
202 # find the offset trend (median offset per axis) in (i, j)-corner
204 xo_median = sorted([xo for xo, _ in offsets])[k]
205 yo_median = sorted([yo for _, yo in offsets])[k]
206 #print('i', i, 'j', j, 'xo_median', xo_median, 'yo_median', yo_median)
208 # choose CORNER_CANDIDATES blocks closest to trend in (i, j)-corner
209 u = numpy.array(offsets, numpy.double)
210 v = numpy.array([xo_median, yo_median], numpy.double)
211 dist = numpy.sum(numpy.square(u - v[numpy.newaxis, :]), 1)
212 corner_candidates.append(
214 [(dist[i],) + blocks[i] for i in range(len(blocks))]
215 )[:CORNER_CANDIDATES]
217 p_all = numpy.stack(p_all, 1)
218 q_all = numpy.stack(q_all, 1)
220 # try all combinations of the corner candidates
221 print('try corner candidates')
222 p = numpy.zeros((2, 4), numpy.double)
223 q = numpy.zeros((2, 4), numpy.double)
226 best_p = None # for diag
227 for _, xf00, yf00, xf10, yf10 in corner_candidates[0]:
232 for _, xf01, yf01, xf11, yf11 in corner_candidates[1]:
237 for _, xf02, yf02, xf12, yf12 in corner_candidates[2]:
242 for _, xf03, yf03, xf13, yf13 in corner_candidates[3]:
248 A = perspective.calc_transform(p, q)
251 q_all - perspective.apply_transform_multi(A, p_all)
254 if best_dist is None or dist < best_dist:
257 best_p = numpy.copy(p) # for diag
260 out_image1 = perspective.remap_image(best_A, image1)
263 xf0, yf0 = numpy.floor(best_p[:, i]).astype(numpy.int32)
266 max(yf0 - 21, 0):max(yf0 + 22, 0),
267 max(xf0 - 1, 0):max(xf0 + 2, 0),
271 max(yf0 - 1, 0):max(yf0 + 2, 0),
272 max(xf0 - 21, 0):max(xf0 + 22, 0),
277 max(yf0 - 21, 0):max(yf0 + 22, 0),
278 max(xf0 - 1, 0):max(xf0 + 2, 0),
282 max(yf0 - 1, 0):max(yf0 + 2, 0),
283 max(xf0 - 21, 0):max(xf0 + 22, 0),
287 sys.stderr.write(f'write {out_jpg0:s}\n')
288 gamma.write_image(out_jpg0, image0)
290 sys.stderr.write(f'write {out_jpg1:s}\n')
291 gamma.write_image(out_jpg1, out_image1)