function [vx, vy] = make_region(convex)
%MAKE_REGION  Random 2D sampling region, returned as CCW vertices (VX, VY).
%   MAKE_REGION(true)  — a convex polygon (vertices on a random ellipse).
%   MAKE_REGION(false) — a non-convex star polygon (star-shaped about the
%                        origin, so the origin is a valid interior start point).
if nargin < 1 || isempty(convex)
    convex = true;
end

if convex
    [vx, vy] = convex_region();
else
    [vx, vy] = star_region();
end

% CCW so the interior is to the left of each edge.
if signed_area(vx, vy) < 0
    vx = vx(end:-1:1);
    vy = vy(end:-1:1);
end
end

function [vx, vy] = convex_region()
% Vertices at increasing angles on an anisotropic, randomly rotated ellipse.
% Points taken in angular order on the ellipse are always in convex position,
% so the polygon is convex by construction — no convhull backend needed (the
% browser worker may not have finished loading it when the figure view
% auto-runs the script). Even angular slots + bounded jitter keep the angles
% ordered and the edges non-degenerate while still random.
m = 6 + randi(4);                       % 7..10 vertices
slot = 2 * pi / m;
ang = (0:m - 1).' * slot + (rand(m, 1) - 0.5) * slot * 0.8;
ex = 1.4 * cos(ang);                    % on an ellipse (anisotropic)
ey = 1.0 * sin(ang);
phi = 2 * pi * rand;                    % random orientation
vx = cos(phi) * ex - sin(phi) * ey;
vy = sin(phi) * ex + cos(phi) * ey;
end

function [vx, vy] = star_region()
% Spikes at sorted angles, alternating outer/inner radius. Small angle jitter
% keeps the angles ordered, so the polygon stays simple and star-shaped about
% the origin (lines through interior points can still leave and re-enter it).
spikes = 4 + randi(4);          % 5..8 spikes
k = 2 * spikes;
step = 2 * pi / k;
ang = (0:k - 1) * step + (rand(1, k) - 0.5) * step * 0.6;
r = zeros(1, k);
r(1:2:k) = 0.9 + 0.4 * rand(1, numel(1:2:k));    % outer
r(2:2:k) = 0.35 + 0.2 * rand(1, numel(2:2:k));   % inner
vx = (1.3 * r .* cos(ang)).';
vy = (1.0 * r .* sin(ang)).';
end

function A = signed_area(vx, vy)
%SIGNED_AREA  Shoelace area; positive when the vertices run counterclockwise.
n = numel(vx);
A = 0;
for i = 1:n
    j = mod(i, n) + 1;
    A = A + (vx(i) * vy(j) - vx(j) * vy(i));
end
A = A / 2;
end