Search code examples
pythonarraysalgorithmmathgeometry

How can I check if all given points (in space) lie on the same line?


I need to implement a function that takes coordinates of any number of points as input data and return True or False depending on whether these points lie on the same line or not.

I use Python to solve this problem and now I have the following implementation:

def are_colinear(points, tolerance): # variable "points" is a list of lists with points coordinates
    for i in range(2, len(points)):
        if (points[i][0] - points[i-2][0])/(points[i-1][0] - points[i-2][0]) - (points[i][1] - points[i-2][1])/(points[i-1][1] - points[i-2][1]) < tolerance and \
           (points[i][1] - points[i-2][1])/(points[i-1][1] - points[i-2][1]) - (points[i][2] - points[i-2][2])/(points[i-1][2] - points[i-2][2]) < tolerance:
            continue
        else:
            return False
    return True

This method is based on the equation of a line passing through two points:

enter image description here

The weakness of this approach is that it raises an error if you want to check points belonging to the same plane (one of three coordinates always equals to zero in this case and one of denominators is zero because of this). I need a better implementation. Thank you in advance!


Solution

  • Take any one of your coordinates, take it to be your new origin, translating all coordinates accordingly. Now, treat each coordinate as a position vector. Normalize each vector. Now, if any two vectors are parallel, their dot product is 1. In fact, they are the same vector. If two vectors are antiparallel, their dot product is -1 and one vector is the negation of the other.

    When you expand it all out, you'll find that you don't need to do any divisions, can avoid any square roots, and don't have special edge cases to handle.

    1 ?= abs(dot(norm(u), norm(v)) 
       = abs(dot(u, v) / (mag(u) * mag(v)))
       = dot(u, v)^2 / (mag(u) * mag(v))^2
    
    1 = (ux*vx + uy*vy + uz*vz)^2 / (sqrt(ux^2 + uv^2 + uz^2) * sqrt(vx^2 + vy^2 + vz^2))^2
    1 = (ux*vx + uy*vy + uz*vz)^2 / ((ux^2 + uv^2 + uz^2) * (vx^2 + vy^2 + vz^2))
    
    (ux^2 + uv^2 + uz^2) * (vx^2 + vy^2 + vz^2) = (ux*vx + uy*vy + uz*vz)^2
    

    This is pretty easy to code up:

    def are_parallel(points):
        points = set(points)  # dedupe
        if len(points) < 3:
            return True
    
        xo, yo, zo = points.pop()  # extract origin
    
        translated_points = [(x-xo, y-yo, z-zo) for x, y, z in points]
        ux, uy, uz = translated_points[0]  # first/reference vector
        u_mag_sqr = ux**2 + uy**2 + uz**2
    
        for vx, vy, vz in translated_points[1:]:  # rest of vectors
            v_mag_sqr = vx**2 + vy**2 + vz**2
            uv_dot_sqr = (ux*vx + uy*vy + uz*vz)**2
            if u_mag_sqr * v_mag_sqr != uv_dot_sqr:
                return False
        return True
    

    Again, its worth emphasizing that this avoids division, square roots, or anything that would introduce floating-point comparison to what could otherwise be integer coordinates, its faster because its just multiplications and additions, and it doesn't have weird edge cases around specific classes of coordinates.