projPointOnCylinder

PROJPOINTONCYLINDER Project a 3D point onto a cylinder.

   PTPROJ = projPointOnCircle3d(PT, CYLINDER).
   Computes the projection of 3D point PT onto the CYLINDER. 
   
   Point PT is a 1-by-3 array, and CYLINDER is a 1-by-7 array.
   Result PTPROJ is a 1-by-3 array, containing the coordinates of the
   projection of PT onto the CYLINDER.

   PTPROJ = projPointOnCircle3d(..., OPT)
   with OPT = 'open' (0) (default) or 'closed' (1), specify if the bases 
   of the cylinder should be included.

   Example
       demoProjPointOnCylinder

   See also
       projPointOnLine3d, projPointOnPlane, projPointOnCircle3d

0001 function ptProj = projPointOnCylinder(pt, cyl, varargin)
0002 %PROJPOINTONCYLINDER Project a 3D point onto a cylinder.
0003 %
0004 %   PTPROJ = projPointOnCircle3d(PT, CYLINDER).
0005 %   Computes the projection of 3D point PT onto the CYLINDER.
0006 %
0007 %   Point PT is a 1-by-3 array, and CYLINDER is a 1-by-7 array.
0008 %   Result PTPROJ is a 1-by-3 array, containing the coordinates of the
0009 %   projection of PT onto the CYLINDER.
0010 %
0011 %   PTPROJ = projPointOnCircle3d(..., OPT)
0012 %   with OPT = 'open' (0) (default) or 'closed' (1), specify if the bases
0013 %   of the cylinder should be included.
0014 %
0015 %   Example
0016 %       demoProjPointOnCylinder
0017 %
0018 %   See also
0019 %       projPointOnLine3d, projPointOnPlane, projPointOnCircle3d
0020 %
0021 
0022 % ------
0023 % Author: oqilipo
0024 % E-mail: N/A
0025 % Created: 2021-04-17, using R2020b
0026 % Copyright 2021-2024
0027 
0028 parser = inputParser;
0029 addRequired(parser, 'pt', @(x) validateattributes(x, {'numeric'},...
0030     {'size',[1 3],'real','finite','nonnan'}));
0031 addRequired(parser, 'cyl', @(x) validateattributes(x, {'numeric'},...
0032     {'size',[1 7],'real','finite','nonnan'}));
0033 capParValidFunc = @(x) (islogical(x) ...
0034     || isequal(x,1) || isequal(x,0) || any(validatestring(x, {'open','closed'})));
0035 addOptional(parser,'cap','open', capParValidFunc);
0036 parse(parser,pt,cyl,varargin{:});
0037 pt = parser.Results.pt;
0038 cyl = parser.Results.cyl;
0039 cap = lower(parser.Results.cap(1));
0040 
0041 % Radius of the cylinder
0042 cylRadius = cyl(7);
0043 % Height of the cylinder
0044 cylBottom = -Inf;
0045 cylHeight = Inf;
0046 if cap == 'c' || cap == 1
0047     cylBottom = 0;
0048     cylHeight = distancePoints3d(cyl(1:3),cyl(4:6));
0049 end
0050 % Create a transformation for the point into a local cylinder coordinate
0051 % system. Align the cylinder axis with the z axis and translate the
0052 % starting point of the cylinder to the origin.
0053 TFM = createRotationVector3d(cyl(4:6)-cyl(1:3), [0 0 1])*createTranslation3d(-cyl(1:3));
0054 % cylTfm = [transformPoint3d(cyl(1:3), TFM) transformPoint3d(cyl(4:6), TFM) cylRadius];
0055 % cylTfm2 = [0 0 0 0 0 cylHeight, cylRadius];
0056 % assert(ismembertol(cylTfm,cylTfm2,'byRows',1,'DataScale',1e1))
0057 
0058 % Transform the point.
0059 ptTfm = transformPoint3d(pt,TFM);
0060 % Convert the transformed point to cylindrical coordinates.
0061 [ptTheta, ptRadius, ptHeight] = cart2cyl(ptTfm);
0062 
0063 if ptRadius <= cylRadius && (ptHeight <= cylBottom || ptHeight >= cylHeight)
0064     % If point is inside the radius of the cylinder but outside its height
0065     if ptHeight <= cylBottom
0066         ptProj_cyl = [ptTheta, ptRadius, 0];
0067     else
0068         ptProj_cyl = [ptTheta, ptRadius, cylHeight];
0069     end
0070 elseif ptRadius > cylRadius && (ptHeight <= cylBottom || ptHeight >= cylHeight)
0071     % If point is outside the cylinder's radius and height
0072     if ptHeight <= cylBottom
0073         ptProj_cyl = [ptTheta, cylRadius, 0];
0074     else
0075         ptProj_cyl = [ptTheta, cylRadius, cylHeight];
0076     end
0077 elseif ptRadius < cylRadius && (ptHeight > cylBottom && ptHeight < cylHeight)
0078     % If point is inside the cylinder's radius and height
0079     deltaRadius = cylRadius - ptRadius;
0080     deltaHeight = cylHeight - ptHeight;
0081     if (deltaRadius < ptHeight && deltaRadius < deltaHeight) || isinf(cylBottom)
0082         % If the distance to the cylinder's surface is smaller than the
0083         % distance to the top and bottom surfaces.
0084         ptProj_cyl = [ptTheta, cylRadius, ptHeight];
0085     else
0086         if ptHeight < deltaHeight
0087             ptProj_cyl = [ptTheta, ptRadius, 0];
0088         else
0089             ptProj_cyl = [ptTheta, ptRadius, cylHeight];
0090         end
0091     end
0092 elseif ptRadius >= cylRadius && (ptHeight > cylBottom && ptHeight < cylHeight)
0093     % If point is outside the radius of the cylinder and inside its height
0094     ptProj_cyl = [ptTheta, cylRadius, ptHeight];
0095 end
0096 
0097 % Convert the projected point back to Cartesian coordinates
0098 ptProj_cart = cyl2cart(ptProj_cyl);
0099 % Transform the projected point back to the global coordinate system
0100 ptProj = transformPoint3d(ptProj_cart,inv(TFM));
0101 
0102 end