Detect collinear points in a 3d points cloud

Lets say we have a random matrix M with N rows and 3 columns (M = randn(N,3)).

I need to detect all the points that are collinear and save the different groups of collinear points

For now here's my idea:

1.I take two points i and i+1 and I browse in the matrix to find the points that are collinear with them

2. If I find collinear points I save their coordinates in a matrix

3. Then I move to the two next set of points i and i+2 and repeat the same idea

It's super slow espeacilly when I use matrix with over 100000 rows. Does anyone have a more optimized approach? Thanks

[Bjorn Gustavsson]

Well, the only thing I can come up with right away is that you should make sure that you only do the necessary tests. That is if you've tested all points for colinearity with point 1 and point 2 you shouldn't do it with point 2 and 1, etc, further if you've found that point l is colinear with points i and j you shouldn't need to do any further tests with points l and i nor with points l and j. There might be additional clever conditions to figure out...

[Sleh Eddine Brika]

Thanks, good point. I could store the tested index combinations in a matrix, each combination is a row, before start comparing a new set I can check if it was tested before

[John D'Errico]

ANY two points are collinear. That is, it is meaningless to ask if TWO points are collinear. 3 points are the minimum.

Therefore, if you want to identify sets of three points as collinear, then you would need to check every combination of 3 points, from a set that is roughly 100000 points. So...

nchoosek(100000,3)
ans =
           166661666700000

I think I counted the digits there. 166 trillion combinations or so?

If your computer can process 1 million sets of three points per second, that leaves

nchoosek(100000,3)/1e6/60/60/24/365
ans =
           5.2848067827245

So around 5 years of CPU time to finish. And, you are surprised this takes a long time? Even if your computer is 10 or 100 times as fast, that is still a BIG problem.

Could you be more efficient? For example, you might decide to loop over over distinct pair of points. There are only

nchoosek(100000,2)
ans =
                4999950000

such pairs. So 5 billion pairs. Determine the line that passes though those points. Now, test to see if any of the other points fall on that line. If you are efficient about this, you can avoid testing a triple when it could have already been seen before. But even so, this is still a big problem, because projecting 100000 points onto a line to see if they fall on the line in question will not be that fast, even if vectorized. You won't be doing this millions of times per second.

You could read what Loren wrote here:

https://blogs.mathworks.com/loren/2008/06/06/collinearity/

But nothing in that blog mentions the problem of how to search through trillions of triples at a time.

[Sleh Eddine Brika]

Thanks for the comment, I never said that I check the collinearity of two points. I just take 2 points and try to find, in the rest of the matrix, what are the points that belong to that line. I will definitely add some more conditions to try simplifying the calculations. Thank you for the recommendations

[Bjorn Gustavsson]

@John: I was thinking in the line of if one have found a set of points that are co-linear with point i and point j then no pair of points from that set have to be used to search for co-linearity with the rest of the ~1e5 points, that's already done. Then if the point-set contains of some number of line-segments with a reasonable number of points on each one would reduce the number of possible lines to test for - possibly by very much.

@Sleh: Maybe you can find efficient methods for this type of analysis if you look for data-analysis-methods for CERN - as far as I understand they do a lot of particle trajectory calculations from point-detections. They also use a lot of computing power to do this.

[Sleh Eddine Brika]

Thanks Bjorn for the suggestion, good idea, I will take a look to their work

Answer 1

[John D'Errico]

I think your problem is in the breath of the search. If you make this a problem of finding sets of three points, then you have trillions of such sets to search over.

So reduce the problem size, not by throwing away points, but by looking for what factor do sets of three collinear points have in common. Not worded well, but you should see the idea once I get going. I want to reduce a doubly or triply nested problem into a single loop.

For example, consider the point cloud

N = 1e5;
xyz = round(rand(N,3)*1000);
xyz = unique(xyz,'rows');

N = size(xyz,1)
N =
       99996

I've rounded them onto an integer lattice, because in my example, I am hoping to find at least a few triplets of points that are truly collinear. In this exmple, I ended up with a set of 99996 distinct points.

Consider the first point in that point cloud. Is it collinear with any other two points? If it is, then there exists one line that passses through point pairs (1,i) and (1,j). So compute the slope of the lines that pass through all such pairs, not as a slope, of course, but as an angle. This does away with the issue of a vertical line, which would have an infinite slope. And, of course, this must be done in 3-dimensions. Again, lets just try it...

A = [atan2(xyz(2:end,2) - xyz(1,2),xyz(2:end,1) - xyz(1,1)),atan2(xyz(2:end,3) - xyz(1,3),xyz(2:end,1) - xyz(1,1))];

So A represents the inclination of each line that passes through point 1, and every other point in the cloud. Essentially, I've computed the slopes dy/dx amd dz/dx for each such line, but I've represented that in an angular form using atan2.

Now, we need to find rows of A that are the same, to within a tolerance. uniquetol does some of the work for us.

[Au,Ia,Ic] = uniquetol(A,1e-12,'byrows','true');
numel(Ia)
ans =
       99942

Now we are getting someplace. There wer 99942 distinct sets of slopes. But that means there were at least a few sets of points that are collinear, and include point #1. We can find them using the Indicator varisble Ic.

colind = find(accumarray(Ic,1) > 1)
colind =
         464
       45424
       77725
       98191
       99941
       99942

So there were 6 distinct lines passing through point #1, such that the line in question also passes through more than one other point. We can easily now go back and find those sets, saving them away.

Having done this for point #1, we need never look at point #1 again, as we have found every line that includes point # 1 and any other two or more points. This means the array A will get progressively smaller as we run along. But it also means the entire process cn be done using just a single loop over the point which may be in common.

In fact, the above scheme will be pretty efficient, since it is just a single loop. Inside the loop create the array of inclinations, then call uniquetol, and use accuarray to find the sets with common inclinations. This will be eminently doable, even for considerably larger point clouds.

[Sleh Eddine Brika]

Thanks John, really clever approach, I will implement it as you proposed. Thank you

Answer 2

[Jeff Miller]

This is going to be a pretty fuzzy suggestion: Maybe you can speed things up by temporarily throwing away precision. The idea is to reduce precision so that you can group points and have fewer to check.

To have a concrete example, imagine that all of your numbers are in the range 0-1. In a first pass, round all the numbers to the nearest 0.1. Now you will have at most 11^3 distinct rows--way less than 100000. Check these distinct rows and keep track of which triples are collinear at this level of precision. This will tell narrow down the triples of the 100000 rows that are possibly collinear, so you only have to check those in a subsequent step. Maybe add another intermediate step where you round to the nearest 0.01.

Well, I'm not sure whether this would work. It would fail if a collinear triple could be rounded into a non-collinear triple. John?