Given two segments in 3D, compute closest point between them using CGAL-CodePudding

Given two segments in 3D in CGAL, I would like to compute the closest points between one segment and another. These two segment may be anywhere in space.

I have looked in CGAL and there is a function which computes distance between both segments (https://doc.cgal.org/latest/Kernel_23/group__squared__distance__grp.html), and with that I guess I could create two spheres and compute the intersection between them, but this seems slow and cumbersome.

Is there something out of the box?

CodePudding user response：

Sorry to say, there are no out-of-the-box way in the CGAL library to see two points, defining the shortest segment between two line segments in 3D. I'm saying "to see" because these two points are actually computed deeply inside the CGAL::squared_distance function, then the distance between them is calculated, and then these two points are discarded.

If you really need these two points you can modify the Distance_3/Segment_3_Segment_3.h file and add a function with an additional argument:

template <class K>
inline
typename K::FT
squared_distance(const Segment_3<K>& seg1,
                 const Segment_3<K>& seg2,
                 Segment_3<K>& res)
{
  // ... modifications
}

Of course you'll need to modify other code in this file as well.

CodePudding user response：

In general, two skew lines define a unique direction orthogonal to both (just take the cross-product of the direction vectors of the lines). Then if you rotate space to make this direction parallel to Z, you reduce the problem to 2D.

Indeed, by Pythagoras, the squared 3D distance is the sum of the squared XY distance and the squared Z distance. The latter is a constant (the two lines are parallel to XY), so it suffices to minimize the XY distance.

You can perform the rotation as a change of basis, creating an orthogonal reference frame with Z as the computed direction; X and Y can be chosen arbitrarily, but the direction of one of the lines is conveniently X and Y is obtained by orthogonality.

The solution of the 2D problem is defined by

the intersection of the segments (if exists),
one endpoint and its orthogonal projection on the other segment (if exists),
two endpoints.

Keep the closest among these point pairs.