geom3d: enhance drawing of multiple geometries at once

matGeom/geom3d/drawAngleBetweenVectors3d.m

@@ -1,4 +1,4 @@
-function varargout = drawAngleBetweenVectors3d(o, v1, v2, r, varargin)
+function varargout = drawAngleBetweenVectors3d(varargin)
 %DRAWANGLEBETWEENVECTORS3D Draw an arc between 2 vectors.
 %
 %   drawAngleBetweenVectors3d(ORIGIN, VECTOR1, VECTOR2, RADIUS) 
@@ -30,41 +30,44 @@
 % Created: 2020-02-02
 % Copyright 2020-2023
 
-% parse axis handle
-hAx = gca;
-if isAxisHandle(o)
-    hAx = o;
-    o = v1;
-    v1 = v2;
-    v2 = r;
-    r = varargin{1};
-    varargin(1) = [];
+% extract handle of axis to draw on
+[hAx, varargin] = parseAxisHandle(varargin{:});
+
+% retrieve inputs
+if length(varargin) < 4
+    error('Requires at least four input arguments');
 end
+o = varargin{1};
+v1 = varargin{2};
+v2 = varargin{3};
+r = varargin{4};
+varargin(1:4) = [];
 
+% parse optional arguments
 p = inputParser;
 p.KeepUnmatched = true;
-logParValidFunc=@(x) (islogical(x) || isequal(x,1) || isequal(x,0));
-addParameter(p,'ConjugateAngle',false,logParValidFunc);
+logParValidFunc = @(x) (islogical(x) || isequal(x,1) || isequal(x,0));
+addParameter(p, 'ConjugateAngle', false, logParValidFunc);
 parse(p, varargin{:});
-conjugate=p.Results.ConjugateAngle;
-drawOptions=p.Unmatched;
+conjugate = p.Results.ConjugateAngle;
+drawOptions = p.Unmatched;
 
 % Normal of the two vectors
-normal=normalizeVector3d(crossProduct3d(v1, v2));
+normal = normalizeVector3d(crossProduct3d(v1, v2));
 % Align normal with the z axis.
 ROT = createRotationVector3d(normal,[0 0 1]);
 % Align first vector with x axis
 ROTv1 = createRotationVector3d(transformVector3d(v1,ROT),[1 0 0]);
 % Get Euler angles of the arc. 
 % The arc is an flat object. Hence, use the 'ZYZ' convention.
-[PHI, THETA, PSI] = rotation3dToEulerAngles((ROTv1*ROT)','ZYZ');
+[PHI, THETA, PSI] = rotation3dToEulerAngles((ROTv1*ROT)', 'ZYZ');
 % Get angle between the vectors
-angle=rad2deg(vectorAngle3d(v1, v2));
+angle = rad2deg(vectorAngle3d(v1, v2));
 % Draw the arc
 if ~conjugate
-    h = drawCircleArc3d(hAx, [o r [THETA PHI PSI] 0 angle],drawOptions);
+    h = drawCircleArc3d(hAx, [o r [THETA PHI PSI] 0 angle], drawOptions);
 else
-    h = drawCircleArc3d(hAx, [o r [THETA PHI PSI] 0 angle-360],drawOptions);
+    h = drawCircleArc3d(hAx, [o r [THETA PHI PSI] 0 angle-360], drawOptions);
 end
 
 % Format output

matGeom/geom3d/drawArrow3d.m

@@ -1,4 +1,4 @@
-function varargout = drawArrow3d(pos, vec, varargin)
+function varargout = drawArrow3d(varargin)
 %DRAWARROW3D Plot a quiver of 3D arrows.
 %
 %   drawArrow3d(pos, vec) 
@@ -43,15 +43,13 @@
 % Created: 2006-09-14, by Shawn Arseneau
 % Copyright 2006-2023
 
-% Check if first argument is an axes handle
-if numel(pos) == 1 && ishghandle(pos, 'axes')
-    hAx = pos;
-    pos=vec;
-    vec=varargin{1};
-    varargin(1)=[];
-else
-    hAx = gca;
-end
+% extract handle of axis to draw on
+[ax, varargin] = parseAxisHandle(varargin{:});
+
+% retrieve positions and vectors
+pos = varargin{1};
+vec = varargin{2};
+varargin(1:2) = [];
 
 numArrows = size(pos,1);
 if numArrows ~= size(vec,1)
@@ -80,14 +78,23 @@
 if numel(stemRatio) == 1; stemRatio = repmat(stemRatio,numArrows,1); end
 arrowRadius = p.Results.arrowRadius;
 if numel(arrowRadius) == 1; arrowRadius = repmat(arrowRadius,numArrows,1); end
-drawOptions=p.Unmatched;
+drawOptions = p.Unmatched;
+
+% save hold state
+holdState = ishold(ax);
+hold(ax, 'on');
+
 
 %% Loop through all arrows and plot in 3D
-hold(hAx,'on')
-qHandle=gobjects(numArrows,1);
-for i=1:numArrows
-    qHandle(i) = drawSingleVector3d(hAx, pos(i,:), vec(i,:), color(i,:), ...
-        stemRatio(i),arrowRadius(i),drawOptions);
+qHandle = gobjects(numArrows,1);
+for i = 1:numArrows
+    qHandle(i) = drawSingleVector3d(ax, pos(i,:), vec(i,:), color(i,:), ...
+        stemRatio(i), arrowRadius(i), drawOptions);
+end
+
+% restore hold state
+if ~holdState
+    hold(ax, 'off');
 end
 
 if nargout > 0
@@ -96,19 +103,19 @@
 
 end
 
-function [valid, color]=validateColor(color,numArrows)
-valid=true;
+function [valid, color] = validateColor(color, numArrows)
+valid = true;
 [arrowRow, arrowCol] = size(color);
-if arrowRow==1
+if arrowRow == 1
     if ischar(color) %in ShortName or LongName color format
-        color=repmat(color,numArrows,1);
+        color = repmat(color,numArrows,1);
     else
-        if arrowCol~=3
+        if arrowCol ~= 3
             error('color in RGBvalue must be of the form 1x3.');
         end
-        color=repmat(color,numArrows,1);
+        color = repmat(color, numArrows, 1);
     end
-elseif arrowRow~=numArrows
+elseif arrowRow ~= numArrows
     error('color in RGBvalue must be of the form Nx3.');
 end
 
@@ -132,16 +139,16 @@
 % Rotate Cylinder
 [row, col] = size(CX); % initial rotation to coincide with x-axis
 
-newEll = transformPoint3d([CX(:), CY(:), CZ(:)],createRotationVector3d([1 0 0],[0 0 -1]));
-CX = reshape(newEll(:,1), row, col);
-CY = reshape(newEll(:,2), row, col);
-CZ = reshape(newEll(:,3), row, col);
+newCyl = transformPoint3d([CX(:), CY(:), CZ(:)], createRotationVector3d([1 0 0],[0 0 -1]));
+CX = reshape(newCyl(:,1), row, col);
+CY = reshape(newCyl(:,2), row, col);
+CZ = reshape(newCyl(:,3), row, col);
 
 [row, col] = size(CX);
-newEll = transformPoint3d([CX(:), CY(:), CZ(:)],createRotationVector3d([1 0 0],vec));
-stemX = reshape(newEll(:,1), row, col);
-stemY = reshape(newEll(:,2), row, col);
-stemZ = reshape(newEll(:,3), row, col);
+newCyl = transformPoint3d([CX(:), CY(:), CZ(:)], createRotationVector3d([1 0 0],vec));
+stemX = reshape(newCyl(:,1), row, col);
+stemY = reshape(newCyl(:,2), row, col);
+stemZ = reshape(newCyl(:,3), row, col);
 
 % Translate cylinder
 stemX = stemX + X;
@@ -159,30 +166,30 @@
 
 % Rotate cone
 [row, col] = size(coneX);
-newEll = transformPoint3d([coneX(:), coneY(:), coneZ(:)],createRotationVector3d([1 0 0],[0 0 -1]));
-coneX = reshape(newEll(:,1), row, col);
-coneY = reshape(newEll(:,2), row, col);
-coneZ = reshape(newEll(:,3), row, col);
+newCyl = transformPoint3d([coneX(:), coneY(:), coneZ(:)], createRotationVector3d([1 0 0], [0 0 -1]));
+coneX = reshape(newCyl(:,1), row, col);
+coneY = reshape(newCyl(:,2), row, col);
+coneZ = reshape(newCyl(:,3), row, col);
 
-newEll = transformPoint3d([coneX(:), coneY(:), coneZ(:)],createRotationVector3d([1 0 0],vec));
-headX = reshape(newEll(:,1), row, col);
-headY = reshape(newEll(:,2), row, col);
-headZ = reshape(newEll(:,3), row, col);
+newCyl = transformPoint3d([coneX(:), coneY(:), coneZ(:)], createRotationVector3d([1 0 0], vec));
+headX = reshape(newCyl(:,1), row, col);
+headY = reshape(newCyl(:,2), row, col);
+headZ = reshape(newCyl(:,3), row, col);
 
 % Translate cone
 % centerline for cylinder: the multiplier is to set the cone 'on the rim' of the cylinder
 V = [0, 0, srho*stemRatio];
-Vp = transformPoint3d(V,createRotationVector3d([1 0 0],[0 0 -1]));
-Vp = transformPoint3d(Vp,createRotationVector3d([1 0 0],vec));
+Vp = transformPoint3d(V, createRotationVector3d([1 0 0], [0 0 -1]));
+Vp = transformPoint3d(Vp, createRotationVector3d([1 0 0], vec));
 headX = headX + Vp(1) + X;
 headY = headY + Vp(2) + Y;
 headZ = headZ + Vp(3) + Z;
 
 % Draw cylinder & cone
 hStem = patch(hAx, surf2patch(stemX, stemY, stemZ), 'FaceColor', color, 'EdgeColor', 'none', drawOptions);
-hold(hAx,'on')
+% hold(hAx,'on')
 hHead = patch(hAx, surf2patch(headX, headY, headZ), 'FaceColor', color, 'EdgeColor', 'none', drawOptions);
 arrowHandle = hggroup(hAx);
-set([hStem, hHead],'Parent',arrowHandle);
+set([hStem, hHead], 'Parent', arrowHandle);
 
 end

matGeom/geom3d/drawBox3d.m

@@ -43,32 +43,41 @@
 zmin = box(:,5);
 zmax = box(:,6);
 
+% save hold state
+holdState = ishold(hAx);
+hold(hAx, 'on');
+
 nBoxes = size(box, 1);
 
-gh=zeros(nBoxes,1);
-for i=1:nBoxes
+gh = zeros(nBoxes,1);
+for i = 1:nBoxes
     % lower face (z=zmin)
-    sh(1) =drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmax(i) ymin(i) zmin(i)], varargin{:});
-    sh(2) =drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmin(i) ymax(i) zmin(i)], varargin{:});
-    sh(3) =drawEdge3d(hAx, [xmax(i) ymin(i) zmin(i)   xmax(i) ymax(i) zmin(i)], varargin{:});
-    sh(4) =drawEdge3d(hAx, [xmin(i) ymax(i) zmin(i)   xmax(i) ymax(i) zmin(i)], varargin{:});
+    sh(1) = drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmax(i) ymin(i) zmin(i)], varargin{:});
+    sh(2) = drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmin(i) ymax(i) zmin(i)], varargin{:});
+    sh(3) = drawEdge3d(hAx, [xmax(i) ymin(i) zmin(i)   xmax(i) ymax(i) zmin(i)], varargin{:});
+    sh(4) = drawEdge3d(hAx, [xmin(i) ymax(i) zmin(i)   xmax(i) ymax(i) zmin(i)], varargin{:});
 
     % front face (y=ymin)
-    sh(5) =drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmin(i) ymin(i) zmax(i)], varargin{:});
-    sh(6) =drawEdge3d(hAx, [xmax(i) ymin(i) zmin(i)   xmax(i) ymin(i) zmax(i)], varargin{:});
-    sh(7) =drawEdge3d(hAx, [xmin(i) ymin(i) zmax(i)   xmax(i) ymin(i) zmax(i)], varargin{:});
+    sh(5) = drawEdge3d(hAx, [xmin(i) ymin(i) zmin(i)   xmin(i) ymin(i) zmax(i)], varargin{:});
+    sh(6) = drawEdge3d(hAx, [xmax(i) ymin(i) zmin(i)   xmax(i) ymin(i) zmax(i)], varargin{:});
+    sh(7) = drawEdge3d(hAx, [xmin(i) ymin(i) zmax(i)   xmax(i) ymin(i) zmax(i)], varargin{:});
 
     % left face (x=xmin)
-    sh(8) =drawEdge3d(hAx, [xmin(i) ymax(i) zmin(i)   xmin(i) ymax(i) zmax(i)], varargin{:});
-    sh(9) =drawEdge3d(hAx, [xmin(i) ymin(i) zmax(i)   xmin(i) ymax(i) zmax(i)], varargin{:});
+    sh(8) = drawEdge3d(hAx, [xmin(i) ymax(i) zmin(i)   xmin(i) ymax(i) zmax(i)], varargin{:});
+    sh(9) = drawEdge3d(hAx, [xmin(i) ymin(i) zmax(i)   xmin(i) ymax(i) zmax(i)], varargin{:});
 
     % the last 3 remaining edges
-    sh(10)=drawEdge3d(hAx, [xmin(i) ymax(i) zmax(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
-    sh(11)=drawEdge3d(hAx, [xmax(i) ymax(i) zmin(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
-    sh(12)=drawEdge3d(hAx, [xmax(i) ymin(i) zmax(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
+    sh(10) = drawEdge3d(hAx, [xmin(i) ymax(i) zmax(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
+    sh(11) = drawEdge3d(hAx, [xmax(i) ymax(i) zmin(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
+    sh(12) = drawEdge3d(hAx, [xmax(i) ymin(i) zmax(i)   xmax(i) ymax(i) zmax(i)], varargin{:});
 
     gh(i) = hggroup(hAx);
-    set(sh,'Parent',gh(i))
+    set(sh, 'Parent', gh(i))
+end
+
+% restore hold state
+if ~holdState
+    hold(hAx, 'off');
 end
 
 if nargout > 0

matGeom/geom3d/drawCapsule.m

@@ -24,7 +24,7 @@
 %   Specifies one or several options using parameter name-value pairs.
 %   Available options are usual drawing options, as well as:
 %   'nPhi'    the number of arcs used for drawing the meridians
-%             (for the semi-spheres and the cylinder(
+%             (for the semi-spheres and the cylinder)
 %   'nTheta'  the number of circles used for drawing the parallels
 %             (only for the semi-spheres at the ends of the capsule)
 %
@@ -85,13 +85,8 @@
 
 %% Input argument processing
 
-% parse axis handle
-if isAxisHandle(varargin{1})
-    hAx = varargin{1};
-    varargin(1) = [];
-else
-    hAx = gca;
-end
+% extract handle of axis to draw on
+[hAx, varargin] = parseAxisHandle(varargin{:});
 
 % input argument representing capsules
 cap = varargin{1};
@@ -152,17 +147,25 @@
     options_cyl(ind:ind+1) = [];
 end
 
-hold on
+% save hold state
+holdState = ishold(hAx);
+hold(hAx, 'on');
+
 if all(cap(1:3) == cap(4:6))
-  % the capsule is only a sphere. take arbitrary axis to be able to plot
-  cap(4:6) = cap(1:3)+eps*([0 0 1]);
-  h1 = 0;
+    % the capsule is only a sphere. take arbitrary axis to be able to plot
+    cap(4:6) = cap(1:3)+eps*([0 0 1]);
+    h1 = 0;
 else
-  h1 = drawCylinder(cap, options_cyl{:});
+    h1 = drawCylinder(cap, options_cyl{:});
 end
 h2 = drawDome(cap([1:3,7]),  (cap(1:3)-cap(4:6)), varargin{:});
 h3 = drawDome(cap([4:6,7]), -(cap(1:3)-cap(4:6)), varargin{:});
 
+% restore hold state
+if ~holdState
+    hold(hAx, 'off');
+end
+
 % return handles
 if nargout == 1
     varargout{1} = [h1, h2, h3];

matGeom/geom3d/drawCircle3d.m

@@ -205,6 +205,11 @@
 nCircles = length(xc);
 h = zeros(nCircles, 1);
 
+% save hold state
+holdState = ishold(hAx);
+hold(hAx, 'on');
+
+% iterate over circles to draw
 for i = 1:nCircles
     % compute position of circle points
     x       = r(i) * cos(t)';
@@ -222,6 +227,10 @@
     h(i) = drawPolyline3d(hAx, circle, options{:});
 end
 
+% restore hold state
+if ~holdState
+    hold(hAx, 'off');
+end
 
 if nargout > 0
     varargout = {h};

matGeom/geom3d/drawCircleArc3d.m

@@ -28,21 +28,35 @@
 varargin(1) = [];
 
 if iscell(arc)
+    % save hold state
+    holdState = ishold(hAx);
+    hold(hAx, 'on');
     h = [];
     for i = 1:length(arc)
         h = [h drawCircleArc3d(hAx, arc{i}, varargin{:})]; %#ok<AGROW>
     end
+    % restore hold state
+    if ~holdState
+        hold(hAx, 'off');
+    end
     if nargout > 0
         varargout = {h};
     end
     return;
 end
 
 if size(arc, 1) > 1
+    % save hold state
+    holdState = ishold(hAx);
+    hold(hAx, 'on');
     h = [];
     for i = 1:size(arc, 1)
         h = [h drawCircleArc3d(hAx, arc(i,:), varargin{:})]; %#ok<AGROW>
     end
+    % restore hold state
+    if ~holdState
+        hold(hAx, 'off');
+    end
     if nargout > 0
         varargout = {h};
     end