Guest User

Untitled

a guest
Jan 30th, 2015
10
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import skimage.measure as skmeasure
  2. import skimage.io as skio
  3. import scipy.stats as stats
  4. import scipy.ndimage as ndimage
  5. import matplotlib.pyplot as plt
  6. import os
  7.  
  8.  
  9. def threshold(im, threshmin=1, threshmax=20, display=False):
  10.  
  11.     image_border = (im == 0) / 2
  12.     water = stats.threshold(im, threshmin, threshmax)
  13.  
  14.     #select only the largest water area, removes all noise and inland lakes and rivers
  15.     labeled_areas, _ = ndimage.label(water)
  16.     mainregion = max(skmeasure.regionprops(labeled_areas), key=lambda r: r.filled_area)
  17.  
  18.     water = labeled_areas == mainregion.label
  19.     water = water + image_border
  20.  
  21.     if display:
  22.         fig, (ax1, ax2) = plt.subplots(1, 2)
  23.  
  24.         ax1.imshow(im, cmap=plt.cm.gray)
  25.         ax1.axis('off')
  26.         ax1.set_title('raw')
  27.  
  28.         ax2.imshow(water, cmap=plt.cm.gray)
  29.         ax2.axis('off')
  30.         ax2.set_title('cleaned')
  31.  
  32.         plt.show()
  33.  
  34.     return water
  35.  
  36.  
  37.  
  38. if __name__ == "__main__":
  39.     #Takes a landsat band 5 image
  40.     filename = os.path.join('testdata', 'LT50400372011313PAC01_B5.TIF')
  41.     im = skio.imread(filename)
  42.     data = threshold(im, display=True)
RAW Paste Data