In [69]:

    

%pylab inline
x = np.linspace(-2*np.pi, 2*np.pi, 100)


    

    




Populating the interactive namespace from numpy and matplotlib






    




WARNING: pylab import has clobbered these variables: ['f']
`%matplotlib` prevents importing * from pylab and numpy

In [70]:

    

f = sin(x)
plot(x, f)


    

    
Out[70]:





[<matplotlib.lines.Line2D at 0x106be9310>]






    

In [116]:

    

# size of structure set
b = 8
# erosion
g = []
for n in range(len(f)):
    start = n-b/2
    stop = n+b/2
    if start < 0:
        continue
    if n+b/2 >= len(f):
        break
    g.append(min(f[start:stop]))
# dilation
h = []
for n in range(len(g)):
    start = n-b/2
    stop = n+b/2
    if start < 0:
        continue
    if n+b/2 >= len(g):
        break
    h.append(max(g[start:stop]))

xg = x[b/2:len(x)-b/2]
xh = xg[b/2:len(xg)-b/2]

figure(figsize=(12,12))
plot(x, f, '-', label='f', linewidth=2)
plot(xg, g, '.', label='erosion', linewidth=2)
plot(xh, h, '*', label='opening', linewidth=2)
legend()
title('Opening')


    

    
Out[116]:





<matplotlib.text.Text at 0x109055fd0>






    

In [104]:

    

def erosion(x, y, structure_size=8):
    """Gray-scale erosion as defined page 668 in Digital Image Processing.
    
    :x: x elements
    :y: y elements
    :structure_size: size of structure, should be even
    :return: tuple of x-array and eroded y-array.
    """
    b = structure_size
    e = []
    for n in range(len(y)):
        start = n-b/2
        stop = n+b/2
        if start < 0:
            continue
        if n+b/2 >= len(y):
            break
        e.append(min(y[start:stop]))
    xe = x[b/2:len(x)-b/2]
    return xe, e


def dilation(x, y, structure_size=8):
    """Gray-scale dilation as defined page 668 in Digital Image Processing.
    
    :x: x elements
    :y: y elements
    :structure_size: size of structure, should be even
    :return: list of tuples with x's and y's
    """
    b = structure_size
    d = []
    for n in range(len(y)):
        start = n-b/2
        stop = n+b/2
        if start < 0:
            continue
        if n+b/2 >= len(y):
            break
        d.append(max(y[start:stop]))
    xd = x[b/2:len(x)-b/2]
    return xd, d


    

In [114]:

    

xd, yd = dilation(x,f)
xe, ye = erosion(xd, yd)
figure(figsize=(12,12))
plot(x,f, '-', label='signal')
plot(xd,yd, '.', label='dilation')
plot(xe,ye, '*', label='closing')
legend()
title('Closing')


    

    
Out[114]:





<matplotlib.text.Text at 0x1065a3050>






    

In [ ]: