In [1]:
%matplotlib inline
import io
import matplotlib.pyplot as plt
import matplotlib.gridspec
from IPython.html.widgets import interact
import skimage.filter
import numpy as np
import scipy.ndimage
import skimage.util 
import skimage.feature
import requests

In [2]:
def plot(frequency, sigma, theta):
    kernel = skimage.filter.gabor_kernel(frequency=frequency, sigma_x=sigma, sigma_y=sigma, theta=theta)
    real, imag = np.real(kernel), np.imag(kernel)
    fig, axes = plt.subplots(1, 2, figsize=(12,6), sharex=True, sharey=True)
    axes[0].imshow(real, cmap='Greys', interpolation='none')
    axes[1].imshow(imag, cmap='Greys', interpolation='none')
interact(plot, frequency=(0.05,0.25, 0.05), sigma=(1,3), theta=(0, np.pi, 0.25*np.pi))


Out[2]:
<function __main__.plot>

In [3]:
url = 'http://argus-public.deltares.nl/sites/jvspeijk/2014/c2/101_Apr.11/1397192405.Fri.Apr.11_05_00_05.GMT.2014.jvspeijk.c2.snap.jpg'
res = requests.get(url)
image = plt.imread(io.BytesIO(res.content), format='jpeg')
image_rgb = plt.imread(open('egmond.jpg'))[12:-12,...]
image = skimage.util.img_as_float(image_rgb).mean(2) #bw float
plt.imshow(image, cmap='Greys_r')


Out[3]:
<matplotlib.image.AxesImage at 0x114513290>

In [4]:
def plot(frequency, sigma, theta, thin):
    imagethin = image[::thin, ::thin]
    kernel = skimage.filter.gabor_kernel(frequency=frequency, sigma_x=sigma, sigma_y=sigma, theta=theta)
    real, imag = skimage.filter.gabor_filter(imagethin, frequency=frequency, sigma_x=sigma, sigma_y=sigma, theta=theta)
    response = np.sqrt(real ** 2 + imag ** 2)
    
    fig = plt.figure(figsize=(13,4.5))
    gs = matplotlib.gridspec.GridSpec(1, 2,width_ratios=[2,3])
    ax0 = plt.subplot(gs[0])
    ax1 = plt.subplot(gs[1])
    # fig, axes = plt.subplots(1, 2, figsize=(13,6))
    ax0.imshow(np.real(kernel), cmap='Greys', interpolation='none')
    ax0.set_title(r'Gabor kernel with frequency %.2f, $\theta$ %.2f' % (sigma, theta))
    ax1.set_title('Response mean: {:.4f} (min: {:.2f}, max: {:.2f})'.format(response.mean(), response.min(), response.max()))
    ax1.imshow(np.sqrt(response), cmap='Greys_r', interpolation='none')
interact(plot, frequency=(0.05,0.25, 0.05), sigma=(1,9), theta=(0, np.pi, 0.05*np.pi), thin=(1,5))



In [5]:
def plot(sigma1):
    thin = 1
    imagethin = image[::thin, ::thin]
    response1 = skimage.filter.gaussian_filter(imagethin, sigma1)
    response2 = skimage.filter.gaussian_filter(imagethin, sigma1 + 1)
    response = response2 - response1
    fig, axes = plt.subplots(1, 1, figsize=(12,6))
    axes.set_title('{:.2f} (min: {:.2f}, max: {:.2f})'.format(response.mean(), response.min(), response.max()))
    im = axes.imshow((np.abs(response)**0.5)*np.sign(response), cmap=matplotlib.cm.BrBG, interpolation='none')
    plt.colorbar(im, ax=axes)
interact(plot, sigma1=(1,21,3))


Out[5]:
<function __main__.plot>

In [6]:
def plot(blocksize, offset):
    plt.imshow(skimage.filter.threshold_adaptive(image,  blocksize, method='gaussian', offset=offset, mode='reflect', param=None), cmap='Greys')
interact(plot, blocksize=(1,9,2), offset=(0,0.01, 0.001))


Out[6]:
<function __main__.plot>

In [7]:
def plot(P, R): 
    plt.imshow(skimage.feature.local_binary_pattern(image, P, R, method='default'), cmap='Set1')
    plt.colorbar()
interact(plot, P=(1,5), R=(0.5, 2.5, 0.5))


Out[7]:
<function __main__.plot>

In [8]:
def plot(sigma): 
    plt.figure(figsize=(13,7))
    plt.imshow(np.log(skimage.feature.corner_shi_tomasi(image, sigma=sigma)), cmap='Greys')
    plt.colorbar()
interact(plot, sigma=(1, 10))


Out[8]:
<function __main__.plot>

In [66]:
def plot(distance, angle):
    props = skimage.measure.regionprops(slic+1, image)
    diss_arr = []
    for prop in props:
        matrix = skimage.feature.greycomatrix((prop.intensity_image * 255).astype('uint8'), [distance], [angle], 256, symmetric=True, normed=True)
        diss = skimage.feature.greycoprops(matrix, 'dissimilarity')[0][0]
        diss_arr.append(diss)
    diss_arr = np.array(diss_arr)
    plt.imshow(diss_arr[slic], cmap='Greens')
interact(plot, distance=(1,10), angle=(0, np.pi, 0.05*np.pi))



In [22]:
import matplotlib.gridspec as gridspec

plt.figure(figsize=(15,5))
gs = gridspec.GridSpec(2, 2,width_ratios=[1,2], height_ratios=[1,3])
ax = plt.subplot(gs[0])
ax.plot([1,3])
ax = plt.subplot(gs[1])
ax.imshow([[1,2]])
ax = plt.subplot(gs[2])
ax.plot([1,3])
ax = plt.subplot(gs[3])
ax.imshow([[1,2]])


Out[22]:
<matplotlib.image.AxesImage at 0x111a80cd0>

In [12]:
matplotlib.figure.Figure.add_subplot?

In [13]:
import skimage.segmentation
slic = skimage.segmentation.slic(image_rgb, n_segments=200, enforce_connectivity=True)
n_pixels = slic.max()

In [57]:
import skimage.feature
import skimage.measure
# trick to fix labels

    
def get_feature(name, channel):
    
    props = skimage.measure.regionprops(slic+1, channel)
    props_arr = np.empty(slic.max()+2, dtype="object")
    for prop in props:
        props_arr[prop.label-1] = prop
    ufunc = np.frompyfunc(lambda x:getattr(x, name, 0), 1, 1)
    return ufunc(props_arr[slic]).astype('double')

prop = skimage.measure.regionprops(slic+1, image)[0]
fig, axes = plt.subplots(3,3, sharex=True, sharey=True, figsize=(20, 15))
for ax, feature in zip(axes[:,0], {'mean_intensity', 'max_intensity', 'min_intensity'}):
    r, g, b = [get_feature(feature, image_rgb[...,i]) for i in range(3)]
    image_slic = np.concatenate([r[:,:,np.newaxis], g[:,:,np.newaxis], b[:,:,np.newaxis]], axis=-1)
    ax.imshow(image_slic/255.0)
    ax.set_title(feature)
    
kwargs_list = [{'cmap': 'YlGnBu'}, {'cmap':'hsv', 'vmin': -(np.pi)/2.0, 'vmax':np.pi/2.0}, {'cmap':'YlGnBu'}]
for ax, feature, kwargs in zip(axes[:,1], 
                           ['solidity', 'orientation', 'eccentricity'], 
                           kwargs_list):
    image_slic = get_feature(feature, image)
    im = ax.imshow(image_slic, **kwargs)
    ax.set_title(feature) 
    plt.colorbar(im, ax=ax)
    ax.set_xlim(0,1400)
    ax.set_ylim(1000, 0)



In [60]:
ax = axes[2,2]
ax.imshow(diss_arr[slic])
fig


Out[60]:

In [68]:
channels = {}
thin = 4
sigma = 4
dsigma = 1
response1 = skimage.filter.gaussian_filter(image, sigma)
response2 = skimage.filter.gaussian_filter(image, sigma + dsigma)
response = response2 - response1
channels['gauss'] = response

real, imag = skimage.filter.gabor_filter(image, frequency=0.15, sigma_x=3, sigma_y=3, theta=1/3.0)
response = np.sqrt(real ** 2 + imag ** 2)
channels['gabor'] = response

kwargs_list = [{'cmap': 'PuOr'}, {'cmap':'Greys'}]
for ax, channel, kwargs in zip(axes[:2,2], 
                           ['gauss', 'gabor'], 
                           kwargs_list):
    image_slic = get_feature('mean_intensity', channels[channel])
    im = ax.imshow(image_slic, **kwargs)
    ax.set_title(channel) 
ax = axes[2,2]    
props = skimage.measure.regionprops(slic+1, image)
diss_arr = []
for prop in props:
    matrix = skimage.feature.greycomatrix((prop.intensity_image * 255).astype('uint8'), [5], [1/3.0], 256, symmetric=True, normed=True)
    diss = skimage.feature.greycoprops(matrix, 'dissimilarity')[0][0]
    diss_arr.append(diss)
diss_arr = np.array(diss_arr)
ax.imshow(diss_arr[slic], cmap='Greens')
ax.set_title('GLCM')
    
fig.tight_layout()
fig


Out[68]:

In [143]:
min([prop.orientation for prop in props])
prop.moments, prop.intensity_image.var()


Out[143]:
(array([[  4.94200000e+03,   1.30963000e+05,   4.05217100e+06,
           1.37216227e+08],
        [  5.02919000e+05,   1.33956330e+07,   4.12886825e+08,
           1.38958525e+10],
        [  6.28528330e+07,   1.67376106e+09,   5.08310003e+10,
           1.67952988e+12],
        [  9.06427692e+09,   2.39845553e+11,   7.13080785e+12,
           2.29740407e+14]]), 4441.2346448353728)

In [117]:
import copy
segmentdata = {}
for name, val in matplotlib.cm.hsv._segmentdata.items():
    val2 = []
    for row in val:
        val2.append((row[0], row[1]/2.0, row[2]/2.0))
    segmentdata[name] = tuple(val2)
    
hsv2 = matplotlib.colors.LinearSegmentedColormap('hsv2', segmentdata)

In [112]:
rgb()


Out[112]:
({'blue': ((0.0, 0.0, 0.0),
   (0.333333, 0.0, 0.0),
   (0.349206, 0.03125, 0.03125),
   (0.507937, 0.5, 0.5),
   (0.84127, 0.5, 0.5),
   (0.857143, 0.46875, 0.46875),
   (1.0, 0.046875, 0.046875)),
  'green': ((0.0, 0.0, 0.0),
   (0.15873, 0.46875, 0.46875),
   (0.174603, 0.5, 0.5),
   (0.507937, 0.5, 0.5),
   (0.666667, 0.03125, 0.03125),
   (0.68254, 0.0, 0.0),
   (1.0, 0.0, 0.0)),
  'red': ((0.0, 0.5, 0.5),
   (0.15873, 0.5, 0.5),
   (0.174603, 0.484375, 0.484375),
   (0.333333, 0.015625, 0.015625),
   (0.349206, 0.0, 0.0),
   (0.666667, 0.0, 0.0),
   (0.68254, 0.015625, 0.015625),
   (0.84127, 0.484375, 0.484375),
   (0.857143, 0.5, 0.5),
   (1.0, 0.5, 0.5))},
 {'blue': ((0.0, 0.0, 0.0),
   (0.333333, 0.0, 0.0),
   (0.349206, 0.015625, 0.015625),
   (0.507937, 0.25, 0.25),
   (0.84127, 0.25, 0.25),
   (0.857143, 0.234375, 0.234375),
   (1.0, 0.0234375, 0.0234375)),
  'green': ((0.0, 0.0, 0.0),
   (0.15873, 0.234375, 0.234375),
   (0.174603, 0.25, 0.25),
   (0.507937, 0.25, 0.25),
   (0.666667, 0.015625, 0.015625),
   (0.68254, 0.0, 0.0),
   (1.0, 0.0, 0.0)),
  'red': ((0.0, 0.25, 0.25),
   (0.15873, 0.25, 0.25),
   (0.174603, 0.2421875, 0.2421875),
   (0.333333, 0.0078125, 0.0078125),
   (0.349206, 0.0, 0.0),
   (0.666667, 0.0, 0.0),
   (0.68254, 0.0078125, 0.0078125),
   (0.84127, 0.2421875, 0.2421875),
   (0.857143, 0.25, 0.25),
   (1.0, 0.25, 0.25))})

In [ ]: