cylinderMesh

CYLINDERMESH Create a 3D mesh representing a cylinder.

   [V, F] = cylinderMesh(CYL)
   Computes vertex coordinates and face vertex indices of a mesh
   representing a 3D cylinder given as [X1 Y1 Z1 X2 Y2 Z2 R].
   
   [V, F] = cylinderMesh(..., OPT)
   with OPT = 'open' (0) (default) or 'closed' (1), specify if the bases 
   of the cylinder should be included.
   
   [V, F] = cylinderMesh(..., NAME, VALUE);
   Specifies one or several options using parameter name-value pairs.
   Available options are:
   'nPerimeter' the number of points represeting the perimeter
   'nRho' the number of circles along the hight

   Example
     % Draw a rotated cylinder
     cyl = [0 0 0 10 20 30 5];
     [v, f] = cylinderMesh(cyl);
     figure;drawMesh(v, f, 'FaceColor', 'r');
     view(3); axis equal;

     % Draw three mutually intersecting cylinders
     p0 = [30 30 30];
     p1 = [90 30 30];
     p2 = [30 90 30];
     p3 = [30 30 90];
     [v1, f1] = cylinderMesh([p0 p1 25]);
     [v2, f2] = cylinderMesh([p0 p2 25]);
     [v3, f3] = cylinderMesh([p0 p3 25],'closed','nPeri',40,'nRho',20);
     figure; hold on;
     drawMesh(v1, f1, 'FaceColor', 'r');
     drawMesh(v2, f2, 'FaceColor', 'g');
     drawMesh(v3, f3, 'FaceColor', 'b');
     view(3); axis equal
     set(gcf, 'renderer', 'opengl')
  
   See also
     drawCylinder, torusMesh, sphereMesh

0001 function varargout = cylinderMesh(cyl, varargin)
0002 %CYLINDERMESH Create a 3D mesh representing a cylinder.
0003 %
0004 %   [V, F] = cylinderMesh(CYL)
0005 %   Computes vertex coordinates and face vertex indices of a mesh
0006 %   representing a 3D cylinder given as [X1 Y1 Z1 X2 Y2 Z2 R].
0007 %
0008 %   [V, F] = cylinderMesh(..., OPT)
0009 %   with OPT = 'open' (0) (default) or 'closed' (1), specify if the bases
0010 %   of the cylinder should be included.
0011 %
0012 %   [V, F] = cylinderMesh(..., NAME, VALUE);
0013 %   Specifies one or several options using parameter name-value pairs.
0014 %   Available options are:
0015 %   'nPerimeter' the number of points represeting the perimeter
0016 %   'nRho' the number of circles along the hight
0017 %
0018 %   Example
0019 %     % Draw a rotated cylinder
0020 %     cyl = [0 0 0 10 20 30 5];
0021 %     [v, f] = cylinderMesh(cyl);
0022 %     figure;drawMesh(v, f, 'FaceColor', 'r');
0023 %     view(3); axis equal;
0024 %
0025 %     % Draw three mutually intersecting cylinders
0026 %     p0 = [30 30 30];
0027 %     p1 = [90 30 30];
0028 %     p2 = [30 90 30];
0029 %     p3 = [30 30 90];
0030 %     [v1, f1] = cylinderMesh([p0 p1 25]);
0031 %     [v2, f2] = cylinderMesh([p0 p2 25]);
0032 %     [v3, f3] = cylinderMesh([p0 p3 25],'closed','nPeri',40,'nRho',20);
0033 %     figure; hold on;
0034 %     drawMesh(v1, f1, 'FaceColor', 'r');
0035 %     drawMesh(v2, f2, 'FaceColor', 'g');
0036 %     drawMesh(v3, f3, 'FaceColor', 'b');
0037 %     view(3); axis equal
0038 %     set(gcf, 'renderer', 'opengl')
0039 %
0040 %   See also
0041 %     drawCylinder, torusMesh, sphereMesh
0042 
0043 % ------
0044 % Author: David Legland
0045 % E-mail: david.legland@inrae.fr
0046 % Created: 2012-10-25, using Matlab 7.9.0.529 (R2009b)
0047 % Copyright 2012-2024 INRA - Cepia Software Platform
0048 
0049 parser = inputParser;
0050 addRequired(parser, 'cyl', @(x) validateattributes(x, {'numeric'},...
0051     {'size',[1 7],'real','finite','nonnan'}));
0052 capParValidFunc = @(x) (islogical(x) ...
0053     || isequal(x,1) || isequal(x,0) || any(validatestring(x, {'open','closed'})));
0054 addOptional(parser,'cap','open', capParValidFunc);
0055 addParameter(parser, 'nPerimeter', 20, @(x) validateattributes(x,{'numeric'},...
0056     {'integer','scalar','>=',4}));
0057 addParameter(parser, 'nRho', 10, @(x) validateattributes(x,{'numeric'},...
0058     {'integer','scalar','>=',2}));
0059 parse(parser,cyl,varargin{:});
0060 cyl=parser.Results.cyl;
0061 cap=lower(parser.Results.cap(1));
0062 NoPP=parser.Results.nPerimeter;
0063 nRho=parser.Results.nRho;
0064 
0065 % extract cylinder data
0066 p1 = cyl(:, 1:3);
0067 p2 = cyl(:, 4:6);
0068 r  = cyl(:, 7);
0069 
0070 % compute length and orientation
0071 [theta, phi, rho] = cart2sph2d(p2 - p1);
0072 
0073 % parametrisation on x
0074 t = linspace(0, 2*pi, NoPP+1);
0075 lx = r * cos(t);
0076 ly = r * sin(t);
0077 
0078 % parametrisation on z
0079 lz = linspace(0, rho, nRho);
0080 
0081 % generate surface grids
0082 x = repmat(lx, [length(lz) 1]);
0083 y = repmat(ly, [length(lz) 1]);
0084 z = repmat(lz', [1 length(t)]);
0085 
0086 % transform points
0087 trans = localToGlobal3d(p1, theta, phi, 0);
0088 [x, y, z] = transformPoint3d(x, y, z, trans);
0089 
0090 % convert to FV mesh
0091 [vertices, faces] = surfToMesh(x, y, z, 'xPeriodic', true);
0092 
0093 % Close cylinder
0094 if cap == 'c' || cap == 1
0095     nR = round(r/(rho/(nRho-1)));
0096     % Base at p1
0097     P1 = circleMesh([0 0 0 r 0 0 0], 'nP',NoPP, 'nR',nR);
0098     P1 = transformMesh(P1, trans);
0099     P1.faces = fliplr(P1.faces);
0100     % Base at p2
0101     P2 = circleMesh([transformPoint3d(p2, inv(trans)) r 0 0 0], 'nP',NoPP, 'nR',nR);
0102     P2 = transformMesh(P2, trans);
0103     % Triangulate the lateral surface of the mesh
0104     latSurf.vertices = vertices;
0105     latSurf.faces = triangulateFaces(faces);
0106     mesh = concatenateMeshes(latSurf, P1, P2);
0107     [vertices, faces] = removeDuplicateVertices(mesh, 1e-8);
0108 end
0109 
0110 % format output
0111 varargout = formatMeshOutput(nargout, vertices, faces);