Add more curvature indices

[otto-link]

def horizontal_curvature(z):
    
    zx, zy = np.gradient(z)
    zxx, zxy = np.gradient(zx)
    _, zyy = np.gradient(zy)
    
    zx2 = zx * zx
    zy2 = zy * zy
    p = zx2 + zy2

    n = zy2 * zxx - 2.0 * zxy * zx * zy + zx2 * zyy
    d = p * np.power(p + 1, 0.5)

    return n / d


def vertical_curvature(z):
    
    zx, zy = np.gradient(z)
    zxx, zxy = np.gradient(zx)
    _, zyy = np.gradient(zy)
    
    zx2 = zx * zx
    zy2 = zy * zy
    p = zx2 + zy2
    
    n = zx2 * zxx + 2.0 * zxy * zx * zy + zy2 * zyy
    d = p * np.power(p + 1, 1.5);

    return n / d


# TODO add pre-filtering
def accumulation_curvature(z):
    kh = horizontal_curvature(z)
    kv = vertical_curvature(z)
    
    return kh * kv