used plt.subplots everywhere

doc/examples/applications/plot_coins_segmentation.py

@@ -16,13 +16,11 @@ from skimage import data
 coins = data.coins()
 hist = np.histogram(coins, bins=np.arange(0, 256))
 
-plt.figure(figsize=(8, 3))
-plt.subplot(121)
-plt.imshow(coins, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.subplot(122)
-plt.plot(hist[1][:-1], hist[0], lw=2)
-plt.title('histogram of grey values')
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 3))
+ax1.imshow(coins, cmap=plt.cm.gray, interpolation='nearest')
+ax1.axis('off')
+ax2.plot(hist[1][:-1], hist[0], lw=2)
+ax2.set_title('histogram of grey values')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -37,17 +35,15 @@ background with the coins:
 
 """
 
-plt.figure(figsize=(6, 3))
-plt.subplot(121)
-plt.imshow(coins > 100, cmap=plt.cm.gray, interpolation='nearest')
-plt.title('coins > 100')
-plt.axis('off')
-plt.subplot(122)
-plt.imshow(coins > 150, cmap=plt.cm.gray, interpolation='nearest')
-plt.title('coins > 150')
-plt.axis('off')
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(6, 3))
+ax1.imshow(coins > 100, cmap=plt.cm.gray, interpolation='nearest')
+ax1.set_title('coins > 100')
+ax1.axis('off')
+ax2.imshow(coins > 150, cmap=plt.cm.gray, interpolation='nearest')
+ax2.set_title('coins > 150')
+ax2.axis('off')
 margins = dict(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0, right=1)
-plt.subplots_adjust(**margins)
+f.subplots_adjust(**margins)
 
 """
 .. image:: PLOT2RST.current_figure
@@ -64,10 +60,10 @@ edge-detector.
 from skimage.filter import canny
 edges = canny(coins/255.)
 
-plt.figure(figsize=(4, 3))
-plt.imshow(edges, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.title('Canny detector')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(edges, cmap=plt.cm.gray, interpolation='nearest')
+ax.axis('off')
+ax.set_title('Canny detector')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -79,10 +75,10 @@ from scipy import ndimage
 
 fill_coins = ndimage.binary_fill_holes(edges)
 
-plt.figure(figsize=(4, 3))
-plt.imshow(fill_coins, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.title('Filling the holes')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(fill_coins, cmap=plt.cm.gray, interpolation='nearest')
+ax.axis('off')
+ax.set_title('Filling the holes')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -93,10 +89,10 @@ objects.
 from skimage import morphology
 coins_cleaned = morphology.remove_small_objects(fill_coins, 21)
 
-plt.figure(figsize=(4, 3))
-plt.imshow(coins_cleaned, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.title('Removing small objects')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(coins_cleaned, cmap=plt.cm.gray, interpolation='nearest')
+ax.axis('off')
+ax.set_title('Removing small objects')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -117,10 +113,10 @@ from skimage.filter import sobel
 
 elevation_map = sobel(coins)
 
-plt.figure(figsize=(4, 3))
-plt.imshow(elevation_map, cmap=plt.cm.jet, interpolation='nearest')
-plt.axis('off')
-plt.title('elevation_map')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(elevation_map, cmap=plt.cm.jet, interpolation='nearest')
+ax.axis('off')
+ax.set_title('elevation_map')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -133,10 +129,10 @@ markers = np.zeros_like(coins)
 markers[coins < 30] = 1
 markers[coins > 150] = 2
 
-plt.figure(figsize=(4, 3))
-plt.imshow(markers, cmap=plt.cm.spectral, interpolation='nearest')
-plt.axis('off')
-plt.title('markers')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(markers, cmap=plt.cm.spectral, interpolation='nearest')
+ax.axis('off')
+ax.set_title('markers')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -147,10 +143,10 @@ starting from the markers determined above:
 """
 segmentation = morphology.watershed(elevation_map, markers)
 
-plt.figure(figsize=(4, 3))
-plt.imshow(segmentation, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.title('segmentation')
+f, ax = plt.subplots(figsize=(4, 3))
+ax.imshow(segmentation, cmap=plt.cm.gray, interpolation='nearest')
+ax.axis('off')
+ax.set_title('segmentation')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -166,16 +162,14 @@ segmentation = ndimage.binary_fill_holes(segmentation - 1)
 labeled_coins, _ = ndimage.label(segmentation)
 image_label_overlay = label2rgb(labeled_coins, image=coins)
 
-plt.figure(figsize=(6, 3))
-plt.subplot(121)
-plt.imshow(coins, cmap=plt.cm.gray, interpolation='nearest')
-plt.contour(segmentation, [0.5], linewidths=1.2, colors='y')
-plt.axis('off')
-plt.subplot(122)
-plt.imshow(image_label_overlay, interpolation='nearest')
-plt.axis('off')
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(6, 3))
+ax1.imshow(coins, cmap=plt.cm.gray, interpolation='nearest')
+ax1.contour(segmentation, [0.5], linewidths=1.2, colors='y')
+ax1.axis('off')
+ax2.imshow(image_label_overlay, interpolation='nearest')
+ax2.axis('off')
 
-plt.subplots_adjust(**margins)
+f.subplots_adjust(**margins)
 
 """
 .. image:: PLOT2RST.current_figure

doc/examples/applications/plot_morphology.py

@@ -32,7 +32,7 @@ from skimage import io
 
 plt.gray()
 phantom = img_as_ubyte(io.imread(data_dir+'/phantom.png', as_grey=True))
-plt.imshow(phantom)
+io.imshow(phantom)
 
 """
 .. image:: PLOT2RST.current_figure

doc/examples/applications/plot_rank_filters.py

@@ -44,13 +44,11 @@ from skimage import data
 noisy_image = img_as_ubyte(data.camera())
 hist = np.histogram(noisy_image, bins=np.arange(0, 256))
 
-plt.figure(figsize=(8, 3))
-plt.subplot(1, 2, 1)
-plt.imshow(noisy_image, interpolation='nearest')
-plt.axis('off')
-plt.subplot(1, 2, 2)
-plt.plot(hist[1][:-1], hist[0], lw=2)
-plt.title('Histogram of grey values')
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 3))
+ax1.imshow(noisy_image, interpolation='nearest', cmap=plt.cm.gray)
+ax1.axis('off')
+ax2.plot(hist[1][:-1], hist[0], lw=2)
+ax2.set_title('Histogram of grey values')
 
 """
 
@@ -77,27 +75,24 @@ noisy_image = img_as_ubyte(data.camera())
 noisy_image[noise > 0.99] = 255
 noisy_image[noise < 0.01] = 0
 
-fig = plt.figure(figsize=(10, 7))
+fig, ax = plt.subplots(2, 2, figsize=(10, 7))
+ax1, ax2, ax3, ax4 = ax.ravel()
 
-plt.subplot(2, 2, 1)
-plt.imshow(noisy_image, vmin=0, vmax=255)
-plt.title('Noisy image')
-plt.axis('off')
+ax1.imshow(noisy_image, vmin=0, vmax=255, cmap=plt.cm.gray)
+ax1.set_title('Noisy image')
+ax1.axis('off')
 
-plt.subplot(2, 2, 2)
-plt.imshow(median(noisy_image, disk(1)), vmin=0, vmax=255)
-plt.title('Median $r=1$')
-plt.axis('off')
+ax2.imshow(median(noisy_image, disk(1)), vmin=0, vmax=255, cmap=plt.cm.gray)
+ax2.set_title('Median $r=1$')
+ax2.axis('off')
 
-plt.subplot(2, 2, 3)
-plt.imshow(median(noisy_image, disk(5)), vmin=0, vmax=255)
-plt.title('Median $r=5$')
-plt.axis('off')
+ax3.imshow(median(noisy_image, disk(5)), vmin=0, vmax=255, cmap=plt.cm.gray)
+ax3.set_title('Median $r=5$')
+ax3.axis('off')
 
-plt.subplot(2, 2, 4)
-plt.imshow(median(noisy_image, disk(20)), vmin=0, vmax=255)
-plt.title('Median $r=20$')
-plt.axis('off')
+ax4.imshow(median(noisy_image, disk(20)), vmin=0, vmax=255, cmap=plt.cm.gray)
+ax4.set_title('Median $r=20$')
+ax4.axis('off')
 
 """
 
@@ -119,19 +114,17 @@ image.
 
 from skimage.filter.rank import mean
 
-fig = plt.figure(figsize=[10, 7])
+fig, (ax1, ax2) = plt.subplots(1, 2, figsize=[10, 7])
 
 loc_mean = mean(noisy_image, disk(10))
 
-plt.subplot(1, 2, 1)
-plt.imshow(noisy_image, vmin=0, vmax=255)
-plt.title('Original')
-plt.axis('off')
+ax1.imshow(noisy_image, vmin=0, vmax=255, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(1, 2, 2)
-plt.imshow(loc_mean, vmin=0, vmax=255)
-plt.title('Local mean $r=10$')
-plt.axis('off')
+ax2.imshow(loc_mean, vmin=0, vmax=255, cmap=plt.cm.gray)
+ax2.set_title('Local mean $r=10$')
+ax2.axis('off')
 
 """
 
@@ -155,25 +148,22 @@ noisy_image = img_as_ubyte(data.camera())
 
 bilat = bilateral_mean(noisy_image.astype(np.uint16), disk(20), s0=10, s1=10)
 
-fig = plt.figure(figsize=[10, 7])
+fig, ax = plt.subplots(2, 2, figsize=(10, 7))
+ax1, ax2, ax3, ax4 = ax.ravel()
 
-plt.subplot(2, 2, 1)
-plt.imshow(noisy_image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+ax1.imshow(noisy_image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(2, 2, 3)
-plt.imshow(bilat, cmap=plt.cm.gray)
-plt.title('Bilateral mean')
-plt.axis('off')
+ax2.imshow(bilat, cmap=plt.cm.gray)
+ax2.set_title('Bilateral mean')
+ax2.axis('off')
 
-plt.subplot(2, 2, 2)
-plt.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax3.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
+ax3.axis('off')
 
-plt.subplot(2, 2, 4)
-plt.imshow(bilat[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax4.imshow(bilat[200:350, 350:450], cmap=plt.cm.gray)
+ax4.axis('off')
 
 """
 
@@ -211,31 +201,26 @@ hist = np.histogram(noisy_image, bins=np.arange(0, 256))
 glob_hist = np.histogram(glob, bins=np.arange(0, 256))
 loc_hist = np.histogram(loc, bins=np.arange(0, 256))
 
-plt.figure(figsize=(10, 10))
+f, ax = plt.subplots(3, 2, figsize=(10, 10))
+ax1, ax2, ax3, ax4, ax5, ax6 = ax.ravel()
 
-plt.subplot(321)
-plt.imshow(noisy_image, interpolation='nearest')
-plt.axis('off')
+ax1.imshow(noisy_image, interpolation='nearest', cmap=plt.cm.gray)
+ax1.axis('off')
 
-plt.subplot(322)
-plt.plot(hist[1][:-1], hist[0], lw=2)
-plt.title('Histogram of gray values')
+ax2.plot(hist[1][:-1], hist[0], lw=2)
+ax2.set_title('Histogram of gray values')
 
-plt.subplot(323)
-plt.imshow(glob, interpolation='nearest')
-plt.axis('off')
+ax3.imshow(glob, interpolation='nearest', cmap=plt.cm.gray)
+ax3.axis('off')
 
-plt.subplot(324)
-plt.plot(glob_hist[1][:-1], glob_hist[0], lw=2)
-plt.title('Histogram of gray values')
+ax4.plot(glob_hist[1][:-1], glob_hist[0], lw=2)
+ax4.set_title('Histogram of gray values')
 
-plt.subplot(325)
-plt.imshow(loc, interpolation='nearest')
-plt.axis('off')
+ax5.imshow(loc, interpolation='nearest', cmap=plt.cm.gray)
+ax5.axis('off')
 
-plt.subplot(326)
-plt.plot(loc_hist[1][:-1], loc_hist[0], lw=2)
-plt.title('Histogram of gray values')
+ax6.plot(loc_hist[1][:-1], loc_hist[0], lw=2)
+ax6.set_title('Histogram of gray values')
 
 """
 
@@ -256,17 +241,15 @@ noisy_image = img_as_ubyte(data.camera())
 
 auto = autolevel(noisy_image.astype(np.uint16), disk(20))
 
-fig = plt.figure(figsize=[10, 7])
+fig, (ax1, ax2) = plt.subplots(1, 2, figsize=[10, 7])
 
-plt.subplot(1, 2, 1)
-plt.imshow(noisy_image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+ax1.imshow(noisy_image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(1, 2, 2)
-plt.imshow(auto, cmap=plt.cm.gray)
-plt.title('Local autolevel')
-plt.axis('off')
+ax2.imshow(auto, cmap=plt.cm.gray)
+ax2.set_title('Local autolevel')
+ax2.axis('off')
 
 """
 
@@ -297,7 +280,7 @@ fig, axes = plt.subplots(nrows=3, figsize=(7, 8))
 ax0, ax1, ax2 = axes
 plt.gray()
 
-ax0.imshow(np.hstack((image, loc_autolevel)))
+ax0.imshow(np.hstack((image, loc_autolevel)), cmap=plt.cm.gray)
 ax0.set_title('Original / auto-level')
 
 ax1.imshow(
@@ -326,24 +309,22 @@ noisy_image = img_as_ubyte(data.camera())
 
 enh = enhance_contrast(noisy_image, disk(5))
 
-fig = plt.figure(figsize=[10, 7])
-plt.subplot(2, 2, 1)
-plt.imshow(noisy_image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+fig, ax = plt.subplots(2, 2, figsize=[10, 7])
+ax1, ax2, ax3, ax4 = ax.ravel()
 
-plt.subplot(2, 2, 3)
-plt.imshow(enh, cmap=plt.cm.gray)
-plt.title('Local morphological contrast enhancement')
-plt.axis('off')
+ax1.imshow(noisy_image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(2, 2, 2)
-plt.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax2.imshow(enh, cmap=plt.cm.gray)
+ax2.set_title('Local morphological contrast enhancement')
+ax2.axis('off')
 
-plt.subplot(2, 2, 4)
-plt.imshow(enh[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax3.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
+ax3.axis('off')
+
+ax4.imshow(enh[200:350, 350:450], cmap=plt.cm.gray)
+ax4.axis('off')
 
 """
 
@@ -360,24 +341,22 @@ noisy_image = img_as_ubyte(data.camera())
 
 penh = enhance_contrast_percentile(noisy_image, disk(5), p0=.1, p1=.9)
 
-fig = plt.figure(figsize=[10, 7])
-plt.subplot(2, 2, 1)
-plt.imshow(noisy_image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+fig, ax = plt.subplots(2, 2, figsize=[10, 7])
+ax1, ax2, ax3, ax4 = ax.ravel()
 
-plt.subplot(2, 2, 3)
-plt.imshow(penh, cmap=plt.cm.gray)
-plt.title('Local percentile morphological\n contrast enhancement')
-plt.axis('off')
+ax1.imshow(noisy_image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(2, 2, 2)
-plt.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax2.imshow(penh, cmap=plt.cm.gray)
+ax2.set_title('Local percentile morphological\n contrast enhancement')
+ax2.axis('off')
 
-plt.subplot(2, 2, 4)
-plt.imshow(penh[200:350, 350:450], cmap=plt.cm.gray)
-plt.axis('off')
+ax3.imshow(noisy_image[200:350, 350:450], cmap=plt.cm.gray)
+ax3.axis('off')
+
+ax4.imshow(penh[200:350, 350:450], cmap=plt.cm.gray)
+ax4.axis('off')
 
 """
 
@@ -460,17 +439,15 @@ m = (np.tile(x, (n, 1)) * np.linspace(0.1, 1, n) * 128 + 128).astype(np.uint8)
 radius = 10
 t = rank.otsu(m, disk(radius))
 
-plt.figure()
+f, (ax1, ax2) = plt.subplots(1, 2)
 
-plt.subplot(1, 2, 1)
-plt.imshow(m)
-plt.title('Original')
-plt.axis('off')
+ax1.imshow(m)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(1, 2, 2)
-plt.imshow(m >= t, interpolation='nearest')
-plt.title('Local Otsu ($r=%d$)' % radius)
-plt.axis('off')
+ax2.imshow(m >= t, interpolation='nearest')
+ax2.set_title('Local Otsu ($r=%d$)' % radius)
+ax2.axis('off')
 
 """
 
@@ -501,27 +478,24 @@ opening = minimum(maximum(noisy_image, disk(5)), disk(5))
 grad = gradient(noisy_image, disk(5))
 
 # display results
-fig = plt.figure(figsize=[10, 7])
+fig, ax = plt.subplots(2, 2, figsize=[10, 7])
+ax1, ax2, ax3, ax4 = ax.ravel()
 
-plt.subplot(2, 2, 1)
-plt.imshow(noisy_image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+ax1.imshow(noisy_image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(2, 2, 2)
-plt.imshow(closing, cmap=plt.cm.gray)
-plt.title('Gray-level closing')
-plt.axis('off')
+ax2.imshow(closing, cmap=plt.cm.gray)
+ax2.set_title('Gray-level closing')
+ax2.axis('off')
 
-plt.subplot(2, 2, 3)
-plt.imshow(opening, cmap=plt.cm.gray)
-plt.title('Gray-level opening')
-plt.axis('off')
+ax3.imshow(opening, cmap=plt.cm.gray)
+ax3.set_title('Gray-level opening')
+ax3.axis('off')
 
-plt.subplot(2, 2, 4)
-plt.imshow(grad, cmap=plt.cm.gray)
-plt.title('Morphological gradient')
-plt.axis('off')
+ax4.imshow(grad, cmap=plt.cm.gray)
+ax4.set_title('Morphological gradient')
+ax4.axis('off')
 
 """
 
@@ -651,12 +625,12 @@ for r in e_range:
 
 rec = np.asarray(rec)
 
-plt.figure()
-plt.title('Performance with respect to element size')
-plt.ylabel('Time (ms)')
-plt.title('Element radius')
-plt.plot(e_range, rec)
-plt.legend(['filter.rank.maximum', 'morphology.dilate'])
+f, ax = plt.subplots()
+ax.set_title('Performance with respect to element size')
+ax.set_ylabel('Time (ms)')
+ax.set_xlabel('Element radius')
+ax.plot(e_range, rec)
+ax.legend(['filter.rank.maximum', 'morphology.dilate'])
 
 """
 
@@ -679,12 +653,12 @@ for s in s_range:
 
 rec = np.asarray(rec)
 
-plt.figure()
-plt.title('Performance with respect to image size')
-plt.ylabel('Time (ms)')
-plt.title('Image size')
-plt.plot(s_range, rec)
-plt.legend(['filter.rank.maximum', 'morphology.dilate'])
+f, ax = plt.subplots()
+ax.set_title('Performance with respect to image size')
+ax.set_ylabel('Time (ms)')
+ax.set_xlabel('Image size')
+ax.plot(s_range, rec)
+ax.legend(['filter.rank.maximum', 'morphology.dilate'])
 
 
 """
@@ -714,13 +688,13 @@ for r in e_range:
 
 rec = np.asarray(rec)
 
-plt.figure()
-plt.title('Performance with respect to element size')
-plt.plot(e_range, rec)
-plt.legend(['filter.rank.median', 'filter.median_filter',
-            'scipy.ndimage.percentile'])
-plt.ylabel('Time (ms)')
-plt.title('Element radius')
+f, ax = plt.subplots()
+ax.set_title('Performance with respect to element size')
+ax.plot(e_range, rec)
+ax.legend(['filter.rank.median', 'filter.median_filter',
+           'scipy.ndimage.percentile'])
+ax.set_ylabel('Time (ms)')
+ax.set_xlabel('Element radius')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -729,10 +703,10 @@ Comparison of outcome of the three methods:
 
 """
 
-plt.figure()
-plt.imshow(np.hstack((rc, rctmf, rndi)))
-plt.title('filter.rank.median vs filtermedian_filter vs scipy.ndimage.percentile')
-plt.axis('off')
+f, ax = plt.subplots()
+ax.imshow(np.hstack((rc, rctmf, rndi)))
+ax.set_title('filter.rank.median vs filtermedian_filter vs scipy.ndimage.percentile')
+ax.axis('off')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -755,13 +729,13 @@ for s in s_range:
 
 rec = np.asarray(rec)
 
-plt.figure()
-plt.title('Performance with respect to image size')
-plt.plot(s_range, rec)
-plt.legend(['filter.rank.median', 'filter.median_filter',
-            'scipy.ndimage.percentile'])
-plt.ylabel('Time (ms)')
-plt.title('Image size')
+f, ax = plt.subplots()
+ax.set_title('Performance with respect to image size')
+ax.plot(s_range, rec)
+ax.legend(['filter.rank.median', 'filter.median_filter',
+           'scipy.ndimage.percentile'])
+ax.set_ylabel('Time (ms)')
+ax.set_xlabel('Image size')
 
 """
 .. image:: PLOT2RST.current_figure

doc/examples/plot_canny.py

@@ -35,24 +35,21 @@ edges1 = filter.canny(im)
 edges2 = filter.canny(im, sigma=3)
 
 # display results
-plt.figure(figsize=(8, 3))
+f, (ax1, ax2, ax3) = plt.subplots(nrows=1, ncols=3, figsize=(8, 3))
 
-plt.subplot(131)
-plt.imshow(im, cmap=plt.cm.jet)
-plt.axis('off')
-plt.title('noisy image', fontsize=20)
+ax1.imshow(im, cmap=plt.cm.jet)
+ax1.axis('off')
+ax1.set_title('noisy image', fontsize=20)
 
-plt.subplot(132)
-plt.imshow(edges1, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('Canny filter, $\sigma=1$', fontsize=20)
+ax2.imshow(edges1, cmap=plt.cm.gray)
+ax2.axis('off')
+ax2.set_title('Canny filter, $\sigma=1$', fontsize=20)
 
-plt.subplot(133)
-plt.imshow(edges2, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('Canny filter, $\sigma=3$', fontsize=20)
+ax3.imshow(edges2, cmap=plt.cm.gray)
+ax3.axis('off')
+ax3.set_title('Canny filter, $\sigma=3$', fontsize=20)
 
-plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.9,
+f.subplots_adjust(wspace=0.02, hspace=0.02, top=0.9,
                     bottom=0.02, left=0.02, right=0.98)
 
 plt.show()

doc/examples/plot_contours.py

@@ -29,12 +29,13 @@ r = np.sin(np.exp((np.sin(x)**3 + np.cos(y)**2)))
 contours = measure.find_contours(r, 0.8)
 
 # Display the image and plot all contours found
-plt.imshow(r, interpolation='nearest')
+f, ax = plt.subplots()
+ax.imshow(r, interpolation='nearest', cmap=plt.cm.gray)
 
 for n, contour in enumerate(contours):
-    plt.plot(contour[:, 1], contour[:, 0], linewidth=2)
+    ax.plot(contour[:, 1], contour[:, 0], linewidth=2)
 
-plt.axis('image')
-plt.xticks([])
-plt.yticks([])
+ax.axis('image')
+ax.set_xticks([])
+ax.set_yticks([])
 plt.show()

doc/examples/plot_convex_hull.py

@@ -40,11 +40,9 @@ image[chull] += 1
 # [ 0.  0.  0.  0.  0.  0.  0.  0.  0.]]
 
 
-fig = plt.subplots(figsize=(10, 6))
-plt.subplot(1, 2, 1)
-plt.title('Original picture')
-plt.imshow(original_image, cmap=plt.cm.gray, interpolation='nearest')
-plt.subplot(1, 2, 2)
-plt.title('Transformed picture')
-plt.imshow(image, cmap=plt.cm.gray, interpolation='nearest')
+fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 6))
+ax1.set_title('Original picture')
+ax1.imshow(original_image, cmap=plt.cm.gray, interpolation='nearest')
+ax2.set_title('Transformed picture')
+ax2.imshow(image, cmap=plt.cm.gray, interpolation='nearest')
 plt.show()

doc/examples/plot_corner.py

@@ -29,9 +29,9 @@ image[230:280, 60:110] = 1
 coords = corner_peaks(corner_harris(image), min_distance=5)
 coords_subpix = corner_subpix(image, coords, window_size=13)
 
-plt.gray()
-plt.imshow(image, interpolation='nearest')
-plt.plot(coords[:, 1], coords[:, 0], '.b', markersize=3)
-plt.plot(coords_subpix[:, 1], coords_subpix[:, 0], '+r', markersize=15)
-plt.axis((0, 350, 350, 0))
+f, ax = plt.subplots()
+ax.imshow(image, interpolation='nearest', cmap=plt.cm.gray)
+ax.plot(coords[:, 1], coords[:, 0], '.b', markersize=3)
+ax.plot(coords_subpix[:, 1], coords_subpix[:, 0], '+r', markersize=15)
+ax.axis((0, 350, 350, 0))
 plt.show()

doc/examples/plot_daisy.py

@@ -20,8 +20,9 @@ img = data.camera()
 descs, descs_img = daisy(img, step=180, radius=58, rings=2, histograms=6,
                          orientations=8, visualize=True)
 
-plt.axis('off')
-plt.imshow(descs_img)
+f, ax = plt.subplots()
+ax.axis('off')
+ax.imshow(descs_img)
 descs_num = descs.shape[0] * descs.shape[1]
-plt.title('%i DAISY descriptors extracted:' % descs_num)
+ax.set_title('%i DAISY descriptors extracted:' % descs_num)
 plt.show()

doc/examples/plot_equalize.py

@@ -92,5 +92,5 @@ ax_cdf.set_ylabel('Fraction of total intensity')
 ax_cdf.set_yticks(np.linspace(0, 1, 5))
 
 # prevent overlap of y-axis labels
-plt.subplots_adjust(wspace=0.4)
+f.subplots_adjust(wspace=0.4)
 plt.show()

doc/examples/plot_hog.py

@@ -90,16 +90,16 @@ image = color.rgb2gray(data.lena())
 fd, hog_image = hog(image, orientations=8, pixels_per_cell=(16, 16),
                     cells_per_block=(1, 1), visualise=True)
 
-plt.figure(figsize=(8, 4))
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4))
 
-plt.subplot(121).set_axis_off()
-plt.imshow(image, cmap=plt.cm.gray)
-plt.title('Input image')
+ax1.axis('off')
+ax1.imshow(image, cmap=plt.cm.gray)
+ax1.set_title('Input image')
 
 # Rescale histogram for better display
 hog_image_rescaled = exposure.rescale_intensity(hog_image, in_range=(0, 0.02))
 
-plt.subplot(122).set_axis_off()
-plt.imshow(hog_image_rescaled, cmap=plt.cm.gray)
-plt.title('Histogram of Oriented Gradients')
+ax2.axis('off')
+ax2.imshow(hog_image_rescaled, cmap=plt.cm.gray)
+ax2.set_title('Histogram of Oriented Gradients')
 plt.show()

doc/examples/plot_holes_and_peaks.py

@@ -22,13 +22,13 @@ image = data.moon()
 image = rescale_intensity(image, in_range=(50, 200))
 
 # convenience function for plotting images
-def imshow(image, **kwargs):
-    plt.figure(figsize=(5, 4))
-    plt.imshow(image, **kwargs)
-    plt.axis('off')
+def imshow(image, title, **kwargs):
+    f, ax = plt.subplots(figsize=(5, 4))
+    ax.imshow(image, **kwargs)
+    ax.axis('off')
+    ax.set_title(title)
 
-imshow(image)
-plt.title('original image')
+imshow(image, 'Original image')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -50,8 +50,7 @@ mask = image
 
 filled = reconstruction(seed, mask, method='erosion')
 
-imshow(filled, vmin=image.min(), vmax=image.max())
-plt.title('after filling holes')
+imshow(filled, 'after filling holes',vmin=image.min(), vmax=image.max())
 
 """
 .. image:: PLOT2RST.current_figure
@@ -62,8 +61,8 @@ isolate the dark regions by subtracting the reconstructed image from the
 original image.
 """
 
-imshow(image - filled)
-plt.title('holes')
+imshow(image - filled, 'holes')
+# plt.title('holes')
 
 """
 .. image:: PLOT2RST.current_figure
@@ -78,8 +77,7 @@ intensity instead of the maximum. The remainder of the process is the same.
 seed = np.copy(image)
 seed[1:-1, 1:-1] = image.min()
 rec = reconstruction(seed, mask, method='dilation')
-imshow(image - rec)
-plt.title('peaks')
+imshow(image - rec, 'peaks')
 plt.show()
 
 """

doc/examples/plot_line_hough_transform.py

@@ -77,30 +77,30 @@ image[idx, idx] = 255
 
 h, theta, d = hough_line(image)
 
-plt.figure(figsize=(8, 4))
+f, ax = plt.subplots(1, 3, figsize=(8, 4))
 
-plt.subplot(131)
-plt.imshow(image, cmap=plt.cm.gray)
-plt.title('Input image')
+ax[0].imshow(image, cmap=plt.cm.gray)
+ax[0].set_title('Input image')
+ax[0].axis('image')
 
-plt.subplot(132)
-plt.imshow(np.log(1 + h),
+ax[1].imshow(np.log(1 + h),
            extent=[np.rad2deg(theta[-1]), np.rad2deg(theta[0]),
                    d[-1], d[0]],
            cmap=plt.cm.gray, aspect=1/1.5)
-plt.title('Hough transform')
-plt.xlabel('Angles (degrees)')
-plt.ylabel('Distance (pixels)')
+ax[1].set_title('Hough transform')
+ax[1].set_xlabel('Angles (degrees)')
+ax[1].set_ylabel('Distance (pixels)')
+ax[1].axis('image')
 
-plt.subplot(133)
-plt.imshow(image, cmap=plt.cm.gray)
+ax[2].imshow(image, cmap=plt.cm.gray)
 rows, cols = image.shape
 for _, angle, dist in zip(*hough_line_peaks(h, theta, d)):
     y0 = (dist - 0 * np.cos(angle)) / np.sin(angle)
     y1 = (dist - cols * np.cos(angle)) / np.sin(angle)
-    plt.plot((0, cols), (y0, y1), '-r')
-plt.axis((0, cols, rows, 0))
-plt.title('Detected lines')
+    ax[2].plot((0, cols), (y0, y1), '-r')
+ax[2].axis((0, cols, rows, 0))
+ax[2].set_title('Detected lines')
+ax[2].axis('image')
 
 # Line finding, using the Probabilistic Hough Transform
 
@@ -108,23 +108,22 @@ image = data.camera()
 edges = canny(image, 2, 1, 25)
 lines = probabilistic_hough_line(edges, threshold=10, line_length=5, line_gap=3)
 
-plt.figure(figsize=(8, 3))
+f2, ax = plt.subplots(1, 3, figsize=(8, 3))
 
-plt.subplot(131)
-plt.imshow(image, cmap=plt.cm.gray)
-plt.title('Input image')
+ax[0].imshow(image, cmap=plt.cm.gray)
+ax[0].set_title('Input image')
+ax[0].axis('image')
 
-plt.subplot(132)
-plt.imshow(edges, cmap=plt.cm.gray)
-plt.title('Canny edges')
+ax[1].imshow(edges, cmap=plt.cm.gray)
+ax[1].set_title('Canny edges')
+ax[1].axis('image')
 
-plt.subplot(133)
-plt.imshow(edges * 0)
+ax[2].imshow(edges * 0)
 
 for line in lines:
     p0, p1 = line
-    plt.plot((p0[0], p1[0]), (p0[1], p1[1]))
+    ax[2].plot((p0[0], p1[0]), (p0[1], p1[1]))
 
-plt.title('Probabilistic Hough')
-plt.axis('image')
+ax[2].set_title('Probabilistic Hough')
+ax[2].axis('image')
 plt.show()

doc/examples/plot_medial_transform.py

@@ -55,14 +55,12 @@ skel, distance = medial_axis(data, return_distance=True)
 # Distance to the background for pixels of the skeleton
 dist_on_skel = distance * skel
 
-plt.figure(figsize=(8, 4))
-plt.subplot(121)
-plt.imshow(data, cmap=plt.cm.gray, interpolation='nearest')
-plt.axis('off')
-plt.subplot(122)
-plt.imshow(dist_on_skel, cmap=plt.cm.spectral, interpolation='nearest')
-plt.contour(data, [0.5], colors='w')
-plt.axis('off')
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4))
+ax1.imshow(data, cmap=plt.cm.gray, interpolation='nearest')
+ax1.axis('off')
+ax2.imshow(dist_on_skel, cmap=plt.cm.spectral, interpolation='nearest')
+ax2.contour(data, [0.5], colors='w')
+ax2.axis('off')
 
-plt.subplots_adjust(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0, right=1)
+f.subplots_adjust(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0, right=1)
 plt.show()

doc/examples/plot_otsu.py

@@ -28,20 +28,17 @@ image = camera()
 thresh = threshold_otsu(image)
 binary = image > thresh
 
-plt.figure(figsize=(8, 2.5))
-plt.subplot(1, 3, 1)
-plt.imshow(image, cmap=plt.cm.gray)
-plt.title('Original')
-plt.axis('off')
+f, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(8, 2.5))
+ax1.imshow(image, cmap=plt.cm.gray)
+ax1.set_title('Original')
+ax1.axis('off')
 
-plt.subplot(1, 3, 2, aspect='equal')
-plt.hist(image)
-plt.title('Histogram')
-plt.axvline(thresh, color='r')
+ax2.hist(image)
+ax2.set_title('Histogram')
+ax2.axvline(thresh, color='r')
 
-plt.subplot(1, 3, 3)
-plt.imshow(binary, cmap=plt.cm.gray)
-plt.title('Thresholded')
-plt.axis('off')
+ax3.imshow(binary, cmap=plt.cm.gray)
+ax3.set_title('Thresholded')
+ax3.axis('off')
 
 plt.show()

doc/examples/plot_peak_local_max.py

@@ -25,25 +25,23 @@ image_max = ndimage.maximum_filter(im, size=20, mode='constant')
 coordinates = peak_local_max(im, min_distance=20)
 
 # display results
-plt.figure(figsize=(8, 3))
-plt.subplot(131)
-plt.imshow(im, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('Original')
+f, ax = plt.subplots(1, 3, figsize=(8, 3))
+ax1, ax2, ax3 = ax.ravel()
+ax1.imshow(im, cmap=plt.cm.gray)
+ax1.axis('off')
+ax1.set_title('Original')
 
-plt.subplot(132)
-plt.imshow(image_max, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('Maximum filter')
+ax2.imshow(image_max, cmap=plt.cm.gray)
+ax2.axis('off')
+ax2.set_title('Maximum filter')
 
-plt.subplot(133)
-plt.imshow(im, cmap=plt.cm.gray)
-plt.autoscale(False)
-plt.plot([p[1] for p in coordinates], [p[0] for p in coordinates], 'r.')
-plt.axis('off')
-plt.title('Peak local max')
+ax3.imshow(im, cmap=plt.cm.gray)
+ax3.autoscale(False)
+ax3.plot([p[1] for p in coordinates], [p[0] for p in coordinates], 'r.')
+ax3.axis('off')
+ax3.set_title('Peak local max')
 
-plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.9,
+f.subplots_adjust(wspace=0.02, hspace=0.02, top=0.9,
                     bottom=0.02, left=0.02, right=0.98)
 
 plt.show()

doc/examples/plot_piecewise_affine.py

@@ -35,7 +35,8 @@ out_rows = image.shape[0] - 1.5 * 50
 out_cols = cols
 out = warp(image, tform, output_shape=(out_rows, out_cols))
 
-plt.imshow(out)
-plt.plot(tform.inverse(src)[:, 0], tform.inverse(src)[:, 1], '.b')
-plt.axis((0, out_cols, out_rows, 0))
+f, ax = plt.subplots()
+ax.imshow(out)
+ax.plot(tform.inverse(src)[:, 0], tform.inverse(src)[:, 1], '.b')
+ax.axis((0, out_cols, out_rows, 0))
 plt.show()

doc/examples/plot_pyramid.py

@@ -30,5 +30,6 @@ for p in pyramid[1:]:
     composite_image[i_row:i_row + n_rows, cols:cols + n_cols] = p
     i_row += n_rows
 
-plt.imshow(composite_image)
+f, ax = plt.subplots()
+ax.imshow(composite_image)
 plt.show()

doc/examples/plot_radon_transform.py

@@ -60,22 +60,20 @@ from skimage.transform import radon, rescale
 image = imread(data_dir + "/phantom.png", as_grey=True)
 image = rescale(image, scale=0.4)
 
-plt.figure(figsize=(8, 4.5))
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4.5))
 
-plt.subplot(121)
-plt.title("Original")
-plt.imshow(image, cmap=plt.cm.Greys_r)
+ax1.set_title("Original")
+ax1.imshow(image, cmap=plt.cm.Greys_r)
 
 theta = np.linspace(0., 180., max(image.shape), endpoint=True)
 sinogram = radon(image, theta=theta, circle=True)
-plt.subplot(122)
-plt.title("Radon transform\n(Sinogram)")
-plt.xlabel("Projection angle (deg)")
-plt.ylabel("Projection position (pixels)")
-plt.imshow(sinogram, cmap=plt.cm.Greys_r,
+ax2.set_title("Radon transform\n(Sinogram)")
+ax2.set_xlabel("Projection angle (deg)")
+ax2.set_ylabel("Projection position (pixels)")
+ax2.imshow(sinogram, cmap=plt.cm.Greys_r,
            extent=(0, 180, 0, sinogram.shape[0]), aspect='auto')
 
-plt.subplots_adjust(hspace=0.4, wspace=0.5)
+f.subplots_adjust(hspace=0.4, wspace=0.5)
 plt.show()
 
 """
@@ -103,13 +101,11 @@ error = reconstruction_fbp - image
 print('FBP rms reconstruction error: %.3g' % np.sqrt(np.mean(error**2)))
 
 imkwargs = dict(vmin=-0.2, vmax=0.2)
-plt.figure(figsize=(8, 4.5))
-plt.subplot(121)
-plt.title("Reconstruction\nFiltered back projection")
-plt.imshow(reconstruction_fbp, cmap=plt.cm.Greys_r)
-plt.subplot(122)
-plt.title("Reconstruction error\nFiltered back projection")
-plt.imshow(reconstruction_fbp - image, cmap=plt.cm.Greys_r, **imkwargs)
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4.5))
+ax1.set_title("Reconstruction\nFiltered back projection")
+ax1.imshow(reconstruction_fbp, cmap=plt.cm.Greys_r)
+ax2.set_title("Reconstruction error\nFiltered back projection")
+ax2.imshow(reconstruction_fbp - image, cmap=plt.cm.Greys_r, **imkwargs)
 plt.show()
 
 """
@@ -156,14 +152,12 @@ error = reconstruction_sart - image
 print('SART (1 iteration) rms reconstruction error: %.3g'
       % np.sqrt(np.mean(error**2)))
 
-plt.figure(figsize=(8, 8.5))
-
-plt.subplot(221)
-plt.title("Reconstruction\nSART")
-plt.imshow(reconstruction_sart, cmap=plt.cm.Greys_r)
-plt.subplot(222)
-plt.title("Reconstruction error\nSART")
-plt.imshow(reconstruction_sart - image, cmap=plt.cm.Greys_r, **imkwargs)
+f, ax = plt.subplots(2, 2, figsize=(8, 8.5))
+ax1, ax2, ax3, ax4 = ax.ravel()
+ax1.set_title("Reconstruction\nSART")
+ax1.imshow(reconstruction_sart, cmap=plt.cm.Greys_r)
+ax2.set_title("Reconstruction error\nSART")
+ax2.imshow(reconstruction_sart - image, cmap=plt.cm.Greys_r, **imkwargs)
 
 # Run a second iteration of SART by supplying the reconstruction
 # from the first iteration as an initial estimate
@@ -173,12 +167,10 @@ error = reconstruction_sart2 - image
 print('SART (2 iterations) rms reconstruction error: %.3g'
       % np.sqrt(np.mean(error**2)))
 
-plt.subplot(223)
-plt.title("Reconstruction\nSART, 2 iterations")
-plt.imshow(reconstruction_sart2, cmap=plt.cm.Greys_r)
-plt.subplot(224)
-plt.title("Reconstruction error\nSART, 2 iterations")
-plt.imshow(reconstruction_sart2 - image, cmap=plt.cm.Greys_r, **imkwargs)
+ax3.set_title("Reconstruction\nSART, 2 iterations")
+ax3.imshow(reconstruction_sart2, cmap=plt.cm.Greys_r)
+ax4.set_title("Reconstruction error\nSART, 2 iterations")
+ax4.imshow(reconstruction_sart2 - image, cmap=plt.cm.Greys_r, **imkwargs)
 plt.show()
 
 """

doc/examples/plot_random_walker_segmentation.py

@@ -58,20 +58,17 @@ markers[data > 1.3] = 2
 labels = random_walker(data, markers, beta=10, mode='bf')
 
 # Plot results
-plt.figure(figsize=(8, 3.2))
-plt.subplot(131)
-plt.imshow(data, cmap='gray', interpolation='nearest')
-plt.axis('off')
-plt.title('Noisy data')
-plt.subplot(132)
-plt.imshow(markers, cmap='hot', interpolation='nearest')
-plt.axis('off')
-plt.title('Markers')
-plt.subplot(133)
-plt.imshow(labels, cmap='gray', interpolation='nearest')
-plt.axis('off')
-plt.title('Segmentation')
+f, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(8, 3.2))
+ax1.imshow(data, cmap='gray', interpolation='nearest')
+ax1.axis('off')
+ax1.set_title('Noisy data')
+ax2.imshow(markers, cmap='hot', interpolation='nearest')
+ax2.axis('off')
+ax2.set_title('Markers')
+ax3.imshow(labels, cmap='gray', interpolation='nearest')
+ax3.axis('off')
+ax3.set_title('Segmentation')
 
-plt.subplots_adjust(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0,
+f.subplots_adjust(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0,
                     right=1)
 plt.show()

doc/examples/plot_ransac.py

@@ -45,11 +45,12 @@ line_x = np.arange(-250, 250)
 line_y = model.predict_y(line_x)
 line_y_robust = model_robust.predict_y(line_x)
 
-plt.plot(data[inliers, 0], data[inliers, 1], '.b', alpha=0.6,
-         label='Inlier data')
-plt.plot(data[outliers, 0], data[outliers, 1], '.r', alpha=0.6,
-         label='Outlier data')
-plt.plot(line_x, line_y, '-k', label='Line model from all data')
-plt.plot(line_x, line_y_robust, '-b', label='Robust line model')
-plt.legend(loc='lower left')
+f, ax = plt.subplots()
+ax.plot(data[inliers, 0], data[inliers, 1], '.b', alpha=0.6,
+        label='Inlier data')
+ax.plot(data[outliers, 0], data[outliers, 1], '.r', alpha=0.6,
+        label='Outlier data')
+ax.plot(line_x, line_y, '-k', label='Line model from all data')
+ax.plot(line_x, line_y_robust, '-b', label='Robust line model')
+ax.legend(loc='lower left')
 plt.show()

doc/examples/plot_regionprops.py

@@ -26,7 +26,8 @@ image = rotate(image, angle=15, order=0)
 label_img = label(image)
 regions = regionprops(label_img)
 
-plt.imshow(image)
+f, ax = plt.subplots()
+ax.imshow(image, cmap=plt.cm.gray)
 
 for props in regions:
     y0, x0 = props.centroid
@@ -36,15 +37,14 @@ for props in regions:
     x2 = x0 - math.sin(orientation) * 0.5 * props.minor_axis_length
     y2 = y0 - math.cos(orientation) * 0.5 * props.minor_axis_length
 
-    plt.plot((x0, x1), (y0, y1), '-r', linewidth=2.5)
-    plt.plot((x0, x2), (y0, y2), '-r', linewidth=2.5)
-    plt.plot(x0, y0, '.g', markersize=15)
+    ax.plot((x0, x1), (y0, y1), '-r', linewidth=2.5)
+    ax.plot((x0, x2), (y0, y2), '-r', linewidth=2.5)
+    ax.plot(x0, y0, '.g', markersize=15)
 
     minr, minc, maxr, maxc = props.bbox
     bx = (minc, maxc, maxc, minc, minc)
     by = (minr, minr, maxr, maxr, minr)
-    plt.plot(bx, by, '-b', linewidth=2.5)
+    ax.plot(bx, by, '-b', linewidth=2.5)
 
-plt.gray()
-plt.axis((0, 600, 600, 0))
+ax.axis((0, 600, 600, 0))
 plt.show()

doc/examples/plot_skeleton.py

@@ -47,19 +47,17 @@ image[circle2] = 0
 skeleton = skeletonize(image)
 
 # display results
-plt.figure(figsize=(8, 4.5))
+f, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4.5))
 
-plt.subplot(121)
-plt.imshow(image, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('original', fontsize=20)
+ax1.imshow(image, cmap=plt.cm.gray)
+ax1.axis('off')
+ax1.set_title('original', fontsize=20)
 
-plt.subplot(122)
-plt.imshow(skeleton, cmap=plt.cm.gray)
-plt.axis('off')
-plt.title('skeleton', fontsize=20)
+ax2.imshow(skeleton, cmap=plt.cm.gray)
+ax2.axis('off')
+ax2.set_title('skeleton', fontsize=20)
 
-plt.subplots_adjust(wspace=0.02, hspace=0.02, top=0.98,
+f.subplots_adjust(wspace=0.02, hspace=0.02, top=0.98,
                     bottom=0.02, left=0.02, right=0.98)
 
 plt.show()