docs: def overlapped_convolution(bin_template, bin_image,
                             tollerance=0.5, splits=(4, 2)):
      """
      As each of these images are hold only binary values, and RFFT2 works on
      float64 greyscale values, we can make the convolution more efficient by
      breaking the image up into :splits: sectons. Each one of these sections
      then has its greyscale value adjusted and then stacked.
  
      We then apply the convolution to this 'stack' of images, and adjust the
      resultant position matches.
      """
      th, tw = bin_template.shape
      ih, iw = bin_image.shape
      hs, ws = splits
      h = ih // hs
      w = iw // ws
      count = numpy.count_nonzero(bin_template)
      assert count > 0
      assert h >= th
      assert w >= tw
  
      yoffset = [(i * h, ((i + 1) * h) + (th - 1)) for i in range(hs)]
      xoffset = [(i * w, ((i + 1) * w) + (tw - 1)) for i in range(ws)]
  
      # image_stacks is Origin (x,y), array, z (height in stack)
      image_stacks = [((x1, y1), bin_image[y1:y2, x1:x2], float((count + 1) ** (num)))
                      for num, (x1, x2, y1, y2) in
                      enumerate((x1, x2, y1, y2) for (x1, x2)
                                in xoffset for (y1, y2) in yoffset)]
  
      pad_h = max(i.shape[0] for _, i, _ in image_stacks)
      pad_w = max(i.shape[1] for _, i, _ in image_stacks)
  
      # rfft metrics must be an even size - why ... maths?
      pad_w += pad_w % 2
      pad_h += pad_h % 2
  
      overlapped_image = sum_2d_images(pad_bin_image_to_shape(i, (pad_h, pad_w))
                                       * num for _, i, num in image_stacks)
      #print "Overlap splits %r, Image Size (%d,%d),
      #Overlapped Size (%d,%d)"  % (splits,iw,ih,pad_w,pad_h)
  
      # Calculate the convolution of the FFT's of the overlapped image & template
      convolution_freqs = (rfft2(overlapped_image) *
                           rfft2(bin_template[::-1, ::-1],
                                 overlapped_image.shape))
  
      # Reverse the FFT to find the result overlapped image
      convolution_image = irfft2(convolution_freqs)
      # At this point, the maximum point in convolution_image should be the
      # bottom right (why?) of the area of greatest match
  
      results = set()
      for (x, y), _, num in image_stacks[::-1]:
          test = convolution_image / num
          filtered = ((test >= (count - tollerance)) &
                     (test <= (count + tollerance)))
          match_points = numpy.transpose(numpy.nonzero(filtered))  # bottom right
          for (fy, fx) in match_points:
              if fx < (tw - 1) or fy < (th - 1):
                  continue
              results.add((x + fx - (tw - 1), y + fy - (th - 1)))
          convolution_image %= num
      return list(results)