docs: def bg_compensate(img, sigma, splinepoints, scale):
      '''Reads file, subtracts background. Returns [compensated image, background].'''
      
      from PIL import Image
      import pylab
      from matplotlib.image import pil_to_array
      from centrosome.filter import canny
      import matplotlib
  
      img = Image.open(img)
      if img.mode=='I;16':
          # 16-bit image
          # deal with the endianness explicitly... I'm not sure
          # why PIL doesn't get this right.
          imgdata = np.fromstring(img.tostring(),np.uint8)
          imgdata.shape=(int(imgdata.shape[0]/2),2)
          imgdata = imgdata.astype(np.uint16)
          hi,lo = (0,1) if img.tag.prefix == 'MM' else (1,0)
          imgdata = imgdata[:,hi]*256 + imgdata[:,lo]
          img_size = list(img.size)
          img_size.reverse()
          new_img = imgdata.reshape(img_size)
          # The magic # for maximum sample value is 281
          if 281 in img.tag:
              img = new_img.astype(np.float32) / img.tag[281][0]
          elif np.max(new_img) < 4096:
              img = new_img.astype(np.float32) / 4095.
          else:
              img = new_img.astype(np.float32) / 65535.
      else:
          img = pil_to_array(img)
      
      pylab.subplot(1,3,1).imshow(img, cmap=matplotlib.cm.Greys_r)
      pylab.show()
      
      if len(img.shape)>2:
          raise ValueError('Image must be grayscale')
  
  ## Create mask that will fix problem when image has black areas outside of well
      edges = canny(img, np.ones(img.shape, bool), 2, .1, .3)
      ci = np.cumsum(edges, 0)
      cj = np.cumsum(edges, 1)
      i,j = np.mgrid[0:img.shape[0], 0:img.shape[1]]
      mask = ci > 0
      mask = mask & (cj > 0)
      mask[1:,:] &= (ci[0:-1,:] < ci[-1,j[0:-1,:]])
      mask[:,1:] &= (cj[:,0:-1] < cj[i[:,0:-1],-1])
      
      import time
      t0 = time.clock()
      bg = backgr(img, mask, MODE_AUTO, sigma, splinepoints=splinepoints, scale=scale)
      print("Executed in %f sec" % (time.clock() - t0))
      bg[~mask] = img[~mask]
  
      pylab.subplot(1,3,2).imshow(img - bg, cmap=matplotlib.cm.Greys_r)
      pylab.subplot(1,3,3).imshow(bg, cmap=matplotlib.cm.Greys_r)
      pylab.show()