import Numeric

def c_Lev(a,b):
    '''Unit cost for insertions, deletions, and substitutions.'''
    assert len(a) <= 1 and len(b) <= 1
    assert not (a == '' and b == '')
    if a == b:
        return 0
    else:
        return 1

def c_indel(a,b):
    '''Unit cost for insertions and deletions.  Substiutions are
    effectively disallowed.'''
    assert len(a) <= 1 and len(b) <= 1
    assert not (a == '' and b == '')
    if a == b:
        return 0
    elif a == '' or b == '':
        return 1
    else:
        return 2

def distance_simple(x, y, c=c_Lev):
    '''Straightforward dynamic programming method for calculating
    weighted edit distance (Levenshtein distance by default).
    Requires quadratic time and space.'''

    m = len(x)
    n = len(y)

    d = Numeric.zeros((m+1,n+1))

    for i in xrange(m+1):
        for j in xrange(n+1):
            if i > 0 and j > 0:
                d[i,j] = min(d[i-1,j]   + c(x[i-1], ''),
                             d[i,  j-1] + c('',     y[j-1]),
                             d[i-1,j-1] + c(x[i-1], y[j-1]))
            elif i > 0:
                d[i,j] = d[i-1,j]   + c(x[i-1], '')
            elif j > 0:
                d[i,j] = d[i,  j-1] + c('',     y[j-1])
    print d
    return d[m,n]

def distance(x, y, c=c_Lev):
    '''Dynamic programming method for calculating weighted edit
    distance (Levenshtein distance by default).  Requires quadratic
    time and linear space.'''

    m = len(x)
    n = len(y)
    if m > n:
        x,y = y,x
        m,n = n,m
    assert len(x) == m and m <= n and n == len(y)

    # allocate two vectors instead of a whole matrix
    d1 = Numeric.zeros(m+1)
    d2 = Numeric.zeros(m+1)

    # initialize the first vector
    for i in xrange(1, m+1):
        d1[i] = d1[i-1] + c(x[i-1], '')

    for j in xrange(n):
        y_j = y[j]
        print d1
        # update the second vector
        d2[0] = d1[0] + c('', y_j)
        for i in xrange(1, m+1):
            d2[i] = min(d2[i-1] + c(x[i-1], ''),
                        d1[i]   + c('', y_j),
                        d1[i-1] + c(x[i-1], y_j))
        d1,d2 = d2,d1

    print d1
    return d1[-1]
