Shortest distance from point to line segment in 3D This function returns the shortest distance from a 3D point P to a line segment defined by two 3D points A and B. It projects P onto the line, then checks if this projected point lies between A and B. If it does, then the distance is calculated using this projected point. 

If the projection is outside the line segment, then the distance is calculated to either point A or B, depending on which is closer.

I wanted to get feedback on if this is a clean solution to the problem.

function d = point_seg_dist(P, A, B)
% Returns distance from a 3D point to a line segment
%
% Inputs
%
%   P : vector
%       Position in space
%
%   A : vector
%       Position of line segment endpoint
%
%   B : vector
%       Position of other line segment endpoint
%
% Outputs
%
%   d : double
%       Distance from point to line segment

AP = P - A; % Vector from A to point
AB = B - A; % Vector from A to B

% Project point onto line
P_line = A + dot(AP, AB) / dot(AB, AB) * AB;

if all(A &lt; P_line) &amp;&amp; all(P_line &lt; B)
    % The point projected onto the line is in between A and B 

    % Projection of point onto segment is the same 
    % as projection of point onto line
    P_seg = P_line; 
else
    % The point projected onto the line is outside of A and B

    if all(P_line &lt;= A)
        % The projected point is closer to A
        P_seg = A;  
    else
        % The projected point is closer to B
        P_seg = B;
    end
end

d = norm(P - P_seg); % Distance to line segment

end
Here is the implementation of the solution provided by Wolfram Alpha, as was suggested in the comments of the OP.

The implemented equation is: 

$$ d=\frac{|(\mathbf{x}_2-\mathbf{x}_1)\times(\mathbf{x}_2-\mathbf{x}_1)|}{|\mathbf{x}_1-\mathbf{x}_0|} $$ 

where

$$x_1=(x_1,y_1,z_1)$$
$$x_2=(x_2,y_2,z_2)$$

are two points on the line and \$x_0\$ is a third point.

function d = point_seg_dist(x0, x1, x2)

    % Make sure we are in 3D
    if length(x0) == 2
        x1(3) = 0;
        x2(3) = 0;
        x0(3) = 0;
    end

    % We want column arrays
    if size(x0,1) &lt; size(x0,2)
         x0 = P';
         x1 = x1';
         x2 = x2';
    end

    d = norm(cross(x2-x1 , x1-x0)) / norm(x2-x1);


Apart from a couple of checks on the input, the calculation is on one line and it does not need projections or checks on the projected point's position.