function h=pdegplot(varargin)
%PDEGPLOT Plot a PDE geometry representation
%
%   PDEGPLOT(G) plots the boundary of the geometry specified by G.
%   G can be any one of the following geometry representations supported
%   by the PDE toolbox:
%      - Decomposed Geometry Matrix, DECSG
%      - Geometry definition defined in a MATLAB-file, PDEGEOM
%      - AnalyticGeometry - an object representing DECSG or PDEGEOM
%      - DiscreteGeometry - discrete geometry object created from STL data
%
%   H=PDEGPLOT(G) returns handles to the plotted graphical objects.
%
%   PDEGPLOT(AX,G,__) plots the geometry to the axes object specified by AX
%   and including any of the name-value pairs below.
%
%   H=PDEGPLOT(G, 'Name',Value,...) can take name-value pairs to control
%   transparency and the display region labels on the plot.
%   Valid name-value pairs include:
%
%   Name                    Value/{Default}
%   ----------------------------------------------------------------
%   VertexLabels            {off} | on. If 'on', vertex labels are plotted.
%
%   EdgeLabels              {off} | on. If 'on', edge labels are plotted.
%
%   FaceLabels              {off} | on. If 'on', face labels are plotted.
%
%   CellLabels              {off} | on. If 'on', cell labels are plotted.
%
%   SubdomainLabels         {off} | on. 2-D only, same as FaceLabels.
%                           Use FaceLabels instead.
%
%   FaceAlpha               {1} | scalar in range [0, 1] 3-D only.
%                           face transparency 1 is fully opaque and
%                           0 is completely transparent.
%
%   Lighting                off | {on}. 3-D only. If 'on', lighting effects
%                           are used.
%
%  Example 1:
%      g = @lshapeg;
%      pdegplot(g, 'EdgeLabels','on')
%
%  Example 2: Import a 3D geometry and plot showing the face labels.
%      pdem = createpde();
%      gm = importGeometry(pdem,'Block.stl')
%      pdegplot(gm, 'FaceLabels','on', 'FaceAlpha',0.5)
%
%   See also PATCH, PDEGEOM, DECSG, CREATEPDE, pde.PDEModel,
%            pde.DiscreteGeometry, pde.AnalyticGeometry

%       Copyright 1994-2024 The MathWorks, Inc.

% check for handle to Axes if passed
[ax,varargin,nargin] = axescheck(varargin{:});
% make sure regular axis, not e.g. PolarAxes
if ~isempty(ax) && ~isa(ax,'matlab.graphics.axis.Axes')
    error(message('pde:pdegplot:BadAxesType'));
end
isUserAx = ~isempty(ax);
if nargin > 0
    g = varargin{1};
    g = convertStringsToChars(g);
end

if nargin > 1
    [varargin{:}] = convertStringsToChars(varargin{:});
end

p = inputParser;
defaultPlotLabels = 'off';

addParameter(p,'subdomainLabels', defaultPlotLabels,@isValidOptionOnOff);
addParameter(p,'cellLabels',defaultPlotLabels,@isValidOptionOnOff);
addParameter(p,'faceLabels',defaultPlotLabels,@isValidOptionOnOff);
addParameter(p,'edgeLabels',defaultPlotLabels,@isValidOptionOnOff);
addParameter(p,'vertexLabels',defaultPlotLabels,@isValidOptionOnOff);
addParameter(p,'faceAlpha',1,@isnumeric);
addParameter(p,'Lighting','on',@isValidOptionOnOff);

% parse the optional arguments
parse(p,varargin{2:end});
plotSubLabels = strcmpi(p.Results.subdomainLabels, 'on');
plotCellLabels = strcmpi(p.Results.cellLabels, 'on');
plotFaceLabels = strcmpi(p.Results.faceLabels, 'on');
plotEdgeLabels = strcmpi(p.Results.edgeLabels, 'on');
plotVertexLabels = strcmpi(p.Results.vertexLabels, 'on');
falpha = p.Results.faceAlpha;
showLighting = strcmpi(p.Results.Lighting, 'on');

if isa(g,'fegeometry')
    g = g.GeomRep;
    if isempty(g)
        error(message('pde:pdegplot:NoGeometryModel'));
    end
end

if isa(g,'pde.EquationModel')
    g = g.Geometry;
    if(isempty(g))
        error(message('pde:pdegplot:NoGeometryModel'));
    end
end

if ~isnumeric(falpha) || ~isscalar(falpha) || ~isfinite(falpha) || falpha < 0 || falpha > 1
    error(message('pde:pdeplot:InvalidAlphaOption'))
end

if isa(g, 'pde.DiscreteGeometry')
    if plotSubLabels
        warning(message('pde:pdegplot:NoSubdomainLabels3D'));
    end
    if g.IsTwoD
        hh = plotTwoDDiscreteGeometry(g, plotFaceLabels, plotEdgeLabels, plotVertexLabels, ax, isUserAx);
    else
        hh = plotThreeDGeometry(g, falpha, plotCellLabels, plotFaceLabels, plotEdgeLabels, plotVertexLabels, ax, isUserAx, showLighting);
    end
else
    plotSubLabels = plotSubLabels | plotFaceLabels;
    if isa(g, 'pde.AnalyticGeometry')
        hh = plotTwoDGeometry(g.geom, plotVertexLabels, plotEdgeLabels, plotSubLabels, ax);
    else
        hh = plotTwoDGeometry(g, plotVertexLabels, plotEdgeLabels, plotSubLabels, ax);
    end
end

if nargout==1
    h=hh;
end

end


function hdl = plotTwoDGeometry(g, plotVertexLabels, plotEdgeLabels, plotSubLabels, ax)
nbs=pdeigeom(g);
d=pdeigeom(g,1:nbs);
if isempty(ax)
    ha = newplot;
else
    ha = ax;
end
if ~ishold(ha)
    cla(ha,'reset');
end
fg = ancestor(ha,'figure');
updateTheme(fg);
% Number of points on a segment
n=50;
X=[];
Y=[];
for i=1:nbs
    s=linspace(d(1,i),d(2,i),n);
    [x1,y1]=pdeigeom(g,i,s);
    X=[X x1 NaN]; %#ok
    Y=[Y y1 NaN]; %#ok
end
if((plotEdgeLabels || plotSubLabels) && ~isnumeric(g))
    dgm = getEFGeomMatrixForFunction(g);
else
    dgm = g;
end
hdl = plot(ha,X,Y);
if(plotEdgeLabels)
    axe = get(hdl, 'Parent');
    pdePlotEdgeLabels( dgm, axe, true );
end
if(plotSubLabels)
    axe = get(hdl, 'Parent');
    meshable = true;
    try
        [p1,~,t1]=initmesh(g,'init','on', 'MesherVersion','R2013a');
    catch
        warning(message('pde:pdegplot:FaceLabelsFailedBadGeom'));
        meshable = false;
    end
    if (meshable)
        numSub = size(find(unique([dgm(6,:) dgm(7,:)])), 2);
        pdePlotSubLabels(p1, t1, numSub, axe, true);
    end
end
if(plotVertexLabels)
    if ~isnumeric(g)
        dgm = getVGeomMatrixForFunction(g);
    else
        dgm = g;
    end
    axe = get(hdl, 'Parent');
    agm = pde.AnalyticGeometry(dgm);
    pdePlotVertexLabels( agm, axe );
    delete(agm);
end
ha = get(hdl,'Parent');
axis(ha,'equal')
axis(ha,'tight');
end


function hdl = plotTwoDDiscreteGeometry(g, plotFaceLabels, plotEdgeLabels, plotVertexLabels, ax, isUserAx)
[Ex, Ey, ~] = g.allDisplayEdges();
if isempty(ax)
    ha = newplot;
else
    ha = ax;
end
if ~ishold(ha)
    cla(ha,'reset');
end
fg = ancestor(ha,'figure');
updateTheme(fg);
if ~isempty(Ex) > 0
    hold(ha,'on');
    hdl = plot(ha,Ex(:),Ey(:), 'Clipping','off');
    hold(ha,'off');
end
if plotFaceLabels
    fxyz = g.faceLabelLocations();
    plotGeomLabels2D(ha, fxyz, 'F');
end
if plotEdgeLabels
    fxyz = g.edgeLabelLocations();
    plotGeomLabels2D(ha, fxyz, 'E');
end
if plotVertexLabels
    fxyz = g.vertexLabelLocations();
    plotGeomLabels2D(ha, fxyz, 'V');
end
axis(ha,'equal');
axis(ha,'tight');
view(ha,0,90);
end

function hdl = plotThreeDGeometry(g, falpha, plotCellLabels, plotFaceLabels, plotEdgeLabels, plotVertexLabels, ax, isUserAx, showLighting)
[FaceFacets, FaceFacetVertices] = g.allDisplayFaces();
[bb, ~] = g.boundingBox();
if isempty(ax)
    ha = newplot;
else
    ha = ax;
end
if ~ishold(ha)
    cla(ha,'reset') % clear axis
end
set(ha,'ClippingStyle','rectangle');
fg = ancestor(ha,'figure');
theme = updateTheme(fg);
hold(ha,'on');
plotAxisTriad(bb,ha);
gainsboro = [hex2dec('DC')/255, hex2dec('DC')/255, hex2dec('DC')/255];
% adjust patch color for dark and light modes
if strcmpi(theme.BaseColorStyle,'dark')
    gainsboro = 6*(ones(1,3)-gainsboro);
end
hdl = patch(ha,'Faces',FaceFacets,'Vertices',FaceFacetVertices,'FaceColor',gainsboro, 'EdgeColor','none', 'FaceAlpha',falpha, 'Clipping','off');
if showLighting && isempty(findobj(ha,"Type","Light"))
    light(ha,'Position',[1,-1,1],'Style','infinite');
    lighting(ha,'flat');
    material(hdl,'dull');
end

[Ex, Ey, Ez] = g.allDisplayEdges();
if ~isempty(Ex) > 0
    hold(ha,'on')
    if strcmpi(theme.BaseColorStyle,'dark')
        hdl(end+1) = plot3(ha,Ex(:),Ey(:),Ez(:),'-w', 'Clipping','off');
    else
        hdl(end+1) = plot3(ha,Ex(:),Ey(:),Ez(:),'-k', 'Clipping','off');
    end
    hold(ha,'off')
end
if plotCellLabels
    fxyz = g.cellLabelLocations();
    plotGeomLabels3D(ha, fxyz, 'C');
end
if plotFaceLabels
    fxyz = g.faceLabelLocations();
    plotGeomLabels3D(ha, fxyz, 'F');
end
if plotEdgeLabels
    fxyz = g.edgeLabelLocations();
    plotGeomLabels3D(ha, fxyz, 'E');
end
if plotVertexLabels
    fxyz = g.vertexLabelLocations();
    plotGeomLabels3D(ha, fxyz, 'V');
end
axis(ha,'equal')
axis(ha,'tight');
view(ha,30,30);
end

function plotGeomLabels2D(ha, xyloc, label)
x = xyloc{1};
y = xyloc{2};
labels = label + string(1:numel(x)); labels = cellstr(labels');
if isempty(labels)
    return
end
hold(ha,'on')
% set label color
fg = ancestor(ha,'figure');
theme = updateTheme(fg);
if strcmpi(theme.BaseColorStyle,'dark')
    labColor = uint8([255 255 255 255]');
else
    labColor = uint8([0 0 0 255]');
end
gt = matlab.graphics.primitive.world.Text('VertexData', single([x y zeros(size(x))])', ...
    'String', labels, 'HorizontalAlignment','center', 'Clipping','on',...
    'Interpreter','none','Parent', ha, 'ColorData', labColor);
gt.Font.Name = 'Helvetica';
gt.Font.Size = 10;
hold(ha,'off')
end

function plotGeomLabels3D(ha, xyzloc, label)
x = xyzloc{1};
y = xyzloc{2};
z = xyzloc{3};
labels = label + string(1:numel(x)); labels = cellstr(labels');
if isempty(labels)
    return
end
hold(ha,'on');
% set label color
fg = ancestor(ha,'figure');
theme = updateTheme(fg);
if strcmpi(theme.BaseColorStyle,'dark')
    labColor = uint8([255 255 255 255]');
else
    labColor = uint8([0 0 0 255]');
end
gt = matlab.graphics.primitive.world.Text('VertexData', single([x y z])', ...
    'String', labels, 'HorizontalAlignment','center', 'Clipping','on',...
    'Interpreter','none','Parent', ha, 'ColorData', labColor);
gt.Font.Name = 'Helvetica';
gt.Font.Size = 10;
hold(ha,'off');
end


% Construct a dummy decomposed geometry matrix for geometry file solely
% for the purpose of displaying edge and face labels.
function dgm=getEFGeomMatrixForFunction(geomFunc)

if(~isa(geomFunc,'function_handle'))
    geomFunc = str2func(geomFunc);
end
ne = geomFunc();
edgeParams=geomFunc(1:ne);
s = (edgeParams(1,:) + edgeParams(2,:))/2;
[x,y] = geomFunc(1:ne, s);
dgm = zeros(10,ne);
dgm(1,:) = 2;
dgm(2,:) = x;
dgm(3,:) = x;
dgm(4,:) = y;
dgm(5,:) = y;
dgm(6,:) = edgeParams(3,:);
dgm(7,:) = edgeParams(4,:);
end

function dgm=getVGeomMatrixForFunction(geomFunc)

if(~isa(geomFunc,'function_handle'))
    geomFunc = str2func(geomFunc);
end
ne = geomFunc();
edgeParams=geomFunc(1:ne);
s = edgeParams(1,:);
[xs,ys] = geomFunc(1:ne, s);
s = edgeParams(2,:);
[xe,ye] = geomFunc(1:ne, s);

dgm = zeros(10,ne);
dgm(1,:) = 2;
dgm(2,:) = xs;
dgm(3,:) = xe;
dgm(4,:) = ys;
dgm(5,:) = ye;
dgm(6,:) = edgeParams(3,:);
dgm(7,:) = edgeParams(4,:);
end

function ok = isValidOptionOnOff(opt)
    try
       validatestring(opt,{'on', 'off'}) ;
    catch
       error(message('pde:pdegplot:GenericMustBeOffOn'));   
    end
    ok = true;
end

function theme = updateTheme(fg)
% Respond to theme changes
if isempty(fg.Theme)
    % only run if figure theme is empty
    matlab.graphics.internal.themes.figureUseDesktopTheme(fg);
end
% if figure has a theme, return it, else return light
theme = matlab.graphics.internal.themes.lightTheme; 
if ~isempty(fg.Theme)
    theme = fg.Theme;
end
end