# Profile Matching in a Point Cloud

A point cloud is generated using uniform random function for (x,y,z). As shown on the following figure, a flat intersecting plane (profile) is being investigated that matches as the best (even if not the exact one) a target profile i.e., given at the left-bottom corner. So the question is:

1- How to find such a match of given target 2D point map through point cloud considering the following notes/conditions?
2- What are then the coordinates/orientations/degree of similarity etc?

Note 1: The profile of interest could be anywhere with any rotation along axes and also could be in different shape e.g., a triangle, rectangle, quadrangle etc depending on its location and orientation. At the following demonstration only a simple rectangle one is shown.

Note 2: A tolerance value could be considered as the distance of the points from the profile. To demonstrate this for the following figure suppose a tolerance of 0.01 times the smallest dimension (~1) so tol=0.01. So if we remove the rest and project all remaining points on the plane of the profile being investigated then we will be able to check its similarity with the target profile.

Note 3: A related topic could be found at Point Pattern Recognition.

• @Developer Off topic but what software are you using to generate those plots? – Spacey Nov 22 '11 at 1:56
• @Mohammad I use Python + MatPlotLib to do my researches and to generate the graphs etc. – Developer Nov 22 '11 at 7:18
• @Developer Fantastic - it is through Python, but what do they mean 'Python shell ala Matlab'? – Spacey Nov 22 '11 at 7:35
• How are the point clouds stored? As a set of coordinates for the center of each point or as a volumetric dataset that has non-zero values in the coordinates around the points? – endolith Nov 22 '11 at 17:23
• @endolith All points have coordinates associated as P:{x,y,z}. They are indeed dimensionless points. With some approximation however they could be discretized to one-pixel dimension as 3D arrays. They may incorporate also other attributes (such as weights etc) over coordinates. – Developer Nov 22 '11 at 23:57

This is always going to require a lot of computation, especially if you want to process as many as 2000 points. I'm sure there are highly-optimized solutions for this kind of pattern-matching already, but you have to figure out what it's called in order to find them.

Since you're talking about a point cloud (sparse data) instead of an image, my cross-correlation method doesn't really apply (and would be even worse computationally). Something like RANSAC probably finds a match quickly, but I don't know much about it.

My attempt at a solution:

Assumptions:

• You want to find the best match, not just a loose or "probably correct" match
• Match will have some small amount of error due to noise in measurement or computation
• Source points are coplanar
• All source points must exist in target (= any unmatched point is a mismatch for the entire profile)

So you should be able to take a lot of shortcuts by disqualifying things and decrease computation time. In short:

1. pick three points from source
2. search through target points, finding sets of 3 points with the same shape
3. when a match of 3 points is found, check all the other points in the plane they define to see if they're a close match
4. if more than one match of all points is found, choose the one with the smallest sum of 3D distances error

More detailed:

pick a point from the source for testing s1 = (x1, y1)
Find nearest point in source s2 = (x2, y2)
d12 = (x1-x2)^2 + (y1-y2)^2
Find second nearest point in source s3 = (x3, y3)
d13 = (x1-x3)^2 + (y1-y3)^2
d23 = (x2-x3)^2 + (y2-y3)^2

for all (x,y,z) test points t1 in target:
# imagine s1 and t1 are coincident
for all other points t2 in target:
if distance from test point > d12:
break out of loop and try another t2 point
if distance ≈ d12:
# imagine source is now rotated so that s1 and s2 are collinear with t1 and t2
for all other points t3 in target:
if distance from t1 > d13 or from t2 > d23:
break and try another t3
if distance from t1 ≈ d13 and from t2 ≈ d23:
# Now you've found matching triangles in source and target
# align source so that s1, s2, s3 are coplanar with t1, t2, t3
project all source points onto this target plane
for all other points in source:
find nearest point in target
measure distance from source point to target point
if it's not within a threshold:
break and try a new t3
else:
sum errors of all matched points for this configuration (defined by t1, t2, t3)


Whichever configuration has the least squared error for all other points is the best match

Since we're working with 3 nearest neighbor test points, matching target points can be simplified by checking if they are within some radius. If searching for a radius of 1 from (0, 0), for instance, we can disqualify (2, 0) based on x1 - x2, without calculating the actual Euclidean distance, to speed it up a bit. This assumes that subtraction is faster than multiplication. There are optimized searches based on more arbitrary fixed radius, too.

function is_closer_than(x1, y1, z1, x2, y2, z2, distance):
if abs(x1 - x2) or abs(y1 - y2) or abs(z1 - z2) > distance:
return False
return (x1 - x2)^2 + (y1 - y2)^2 + (z1 - z2)^2 > distance^2 # sqrt is slow


3D Euclidean distance is $d = \sqrt{(x_1 - x_2)^2 + (y_1 - y_2)^2+(z_1 - z_2)^2}$ but square roots are slow. Since we're just comparing relative distance, I think we could drop the sqrt and consider "distance" to be sums of squares.

The minimum calculation time would be if no 2-point matches are found. If there are 2000 points in the target, this would be 2000 * 2000 distance calculations, though many would be disqualified by a subtraction, and the results of previous calculations could be stored so you only have to do $2000 \choose 2$ = 1,999,000.

Actually, since you will need to calculate all of these anyway, whether you find matches or not, and since you only care about nearest neighbors for this step, if you have the memory it's probably better to pre-compute these values using an optimized algorithm. Something like a Delaunay or Pitteway triangulation, where every point in the target is connected to its nearest neighbors. Store those in a table, then look them up for each point when trying to fit the source triangle to one of the target triangles.

There's a lot of calculations involved, but it should be relatively fast since it's only operating on the data, which is sparse, instead of multiplying lots of meaningless zeros together like cross-correlation of volumetric data would involve. This same idea would work for the 2D case if you found the centers of the dots first and stored them as a set of coordinates.

• The first part in your answer is actually a brute-force method seeking nearby points (counting with respect to a threshold) around all possible planes through the point cloud. It is extremely computation intensive e.g., for only 2000 points a number of 2,662,668,000,000 (Formula) distance computation will be required! – Developer Nov 27 '11 at 11:00
• @Developer: Yes, it's going to take a lot of calculations, especially if you have thousands of points. Yes, for 2000 points, if you don't find any planes, you'd end up doing 2,658,673,998,000 calculations. Presumably you would find planes, though, which would reduce the time because it stops as soon as it's found enough points. But anyway, I was thinking about this and probably have a better idea, and I'll change the answer. – endolith Nov 27 '11 at 22:08
• You absolutely got the point completely right. Just to add that the stopping criteria cannot apply even after finding a suitable plane while it might be that there is a much better match, so all possible planes need to be checked. I have already implemented this idea and found even with the help of Fortran for numbers higher than 500 points it is going to be impossible to have experiences with PC. – Developer Nov 27 '11 at 23:00

I would add @mirror2image description on the alternative solution beside RANSAC, you may consider ICP algorithm (iterative closest point), a description can be found here!

I think the next challenge in using this ICP is to define your own cost function, and the starting pose of the target plane with respect to the 3d cloud point data. Some practical approach is to introduce some random noise in the data during the iteration to avoid converging to the false minima. This is the heuristic part I guess you need to design.

Update:

The steps in simplified form is:

1. Find the closest point for each input point.
2. Compute the transform from the input to the target, and then move the input points using the transform.
3. Compute the similarity function (e.g. the distance for each input point w.r.t to its corresponding pair target point).
4. Check the stopping condition.

Iterate the step 1-4.

There is available library that you can consider here! (I have not tried it yet), there is one section on the registration part (including other methods).

• Thanks for the link and suggestion. Such useful ideas help always us as novice researchers to learn thing quicker. I always appreciate more explanations. – Developer Nov 28 '11 at 11:16