+def bandpass(image):
+ _, _, cs = image.shape
+ return numpy.stack(
+ [
+ scipy.ndimage.gaussian_filter(image[:, :, i], CUTOFF1, mode = 'mirror') -
+ scipy.ndimage.gaussian_filter(image[:, :, i], CUTOFF0, mode = 'mirror')
+ for i in range(cs)
+ ],
+ 2
+ )
+
+def correlate(image0_bp, x0, y0, image1_bp, x1, y1):
+ block0 = image0_bp[
+ y0:y0 + BLOCK_SIZE * BLOCKS,
+ x0:x0 + BLOCK_SIZE * BLOCKS,
+ :
+ ]
+ block1 = image1_bp[
+ y1:y1 + BLOCK_SIZE * (BLOCKS + 1),
+ x1:x1 + BLOCK_SIZE * (BLOCKS + 1),
+ :
+ ]
+
+ # note: swapping block1, block0 flips the output (subtracts x and y
+ # from BLOCK_SIZE) and we need this to find matching part of image1
+ corr = numpy.sum(
+ numpy.stack(
+ [
+ scipy.signal.correlate(
+ block1[:, :, i],
+ block0[:, :, i],
+ mode = 'valid'
+ )
+ for i in range(cs)
+ ],
+ 0
+ ),
+ 0
+ )
+ #temp = corr - numpy.mean(corr)
+ #temp /= 10. * numpy.sqrt(numpy.mean(numpy.square(temp)))
+ #gamma.write_image(f'corr_{i:d}_{j:d}.jpg', temp + .5)
+
+ y, x = numpy.unravel_index(numpy.argmax(corr), corr.shape)
+ return (
+ (x, y)
+ if (
+ x >= CUTOFF1 and
+ x <= BLOCK_SIZE - CUTOFF1 and
+ y >= CUTOFF1 and
+ y <= BLOCK_SIZE - CUTOFF1
+ ) else
+ None
+ )
+