docs: def get_adaptive_threshold(threshold_method, image, threshold,
                             mask = None,
                             adaptive_window_size = 10,
                             **kwargs):
      
      """Given a global threshold, compute a threshold per pixel
      
      Break the image into blocks, computing the threshold per block.
      Afterwards, constrain the block threshold to .7 T < t < 1.5 T.
      
      Block sizes must be at least 50x50. Images > 500 x 500 get 10x10
      blocks.
      """
      # for the X and Y direction, find the # of blocks, given the
      # size constraints
      image_size = np.array(image.shape[:2],dtype=int)
      nblocks = image_size // adaptive_window_size
      #
      # Use a floating point block size to apportion the roundoff
      # roughly equally to each block
      #
      increment = ( np.array(image_size,dtype=float) / 
                    np.array(nblocks,dtype=float))
      #
      # Put the answer here
      #
      thresh_out = np.zeros(image_size, image.dtype)
      #
      # Loop once per block, computing the "global" threshold within the
      # block.
      #
      block_threshold = np.zeros([nblocks[0],nblocks[1]])
      for i in range(nblocks[0]):
          i0 = int(i*increment[0])
          i1 = int((i+1)*increment[0])
          for j in range(nblocks[1]):
              j0 = int(j*increment[1])
              j1 = int((j+1)*increment[1])
              block = image[i0:i1,j0:j1]
              block_mask = None if mask is None else mask[i0:i1,j0:j1]
              block_threshold[i,j] = get_global_threshold(
                  threshold_method, 
                  block, mask = block_mask,
                  **kwargs)
      #
      # Use a cubic spline to blend the thresholds across the image to avoid image artifacts
      #
      spline_order = min(3, np.min(nblocks) - 1)
      xStart = int(increment[0] / 2)
      xEnd = int((nblocks[0] - 0.5) * increment[0])
      yStart = int(increment[1] / 2)
      yEnd = int((nblocks[1] - 0.5) * increment[1])
      xtStart = .5
      xtEnd = image.shape[0] - .5
      ytStart = .5
      ytEnd = image.shape[1] - .5
      block_x_coords = np.linspace(xStart,xEnd, nblocks[0])
      block_y_coords = np.linspace(yStart,yEnd, nblocks[1])
      adaptive_interpolation = scipy.interpolate.RectBivariateSpline(
          block_x_coords, block_y_coords, block_threshold,
          bbox = (xtStart, xtEnd, ytStart, ytEnd),
          kx = spline_order, ky = spline_order)
      thresh_out_x_coords = np.linspace(.5, int(nblocks[0] * increment[0]) - .5, thresh_out.shape[0])
      thresh_out_y_coords = np.linspace(.5, int(nblocks[1] * increment[1]) - .5 , thresh_out.shape[1])
  
      thresh_out = adaptive_interpolation(thresh_out_x_coords, thresh_out_y_coords)
      
      return thresh_out