def bezier_point(cps, t):
    """
    Cubic Bezier curve interpolation
    B(t) = (1-t)^3 * P0 + 3t(1-t)^2 * P1 + 3t^2(1-t) * P2 + t^3 * P3
    """
    p = (1 - t) ** 3 * cps[(0), :]
    p += 3 * t * (1 - t) ** 2 * cps[(1), :]
    p += 3 * t ** 2 * (1 - t) * cps[(2), :]
    p += t ** 3 * cps[(3), :]
    return p
