def plane_fit(points):
    """
    This fits an n-dimensional plane to a set of points. See
    http://stackoverflow.com/questions/12299540/plane-fitting-to-4-or-more-xyz-points

    :parameter points:
        An instance of :class:~numpy.ndarray. The number of columns must be
        equal to three.
    :return:
         A point on the plane and the normal to the plane.
    """
    points = numpy.transpose(points)
    points = numpy.reshape(points, (numpy.shape(points)[0], -1))
    assert points.shape[0] < points.shape[1], points.shape
    ctr = points.mean(axis=1)
    x = points - ctr[:, None]
    M = numpy.dot(x, x.T)
    return ctr, numpy.linalg.svd(M)[0][:, -1]