From a55798c599e32cf47db3b1e10739f7cc378e124a Mon Sep 17 00:00:00 2001 From: Kevin D Date: Thu, 8 Oct 2026 23:58:09 -0700 Subject: [PATCH] centralized thesis figure plotting scripts --- test/test_rfSensor.m | 25 ----- test/test_sigmoidSensor.m | 24 ---- test_thesis_plots.m | 230 ++++++++++++++++++++++++++++++++++++++ util/exportThesisFig.m | 53 +++++++-- 4 files changed, 272 insertions(+), 60 deletions(-) create mode 100644 test_thesis_plots.m diff --git a/test/test_rfSensor.m b/test/test_rfSensor.m index a03146d..e6a1eae 100644 --- a/test/test_rfSensor.m +++ b/test/test_rfSensor.m @@ -172,30 +172,5 @@ classdef test_rfSensor < matlab.unittest.TestCase sensor2.plot(pos2(3), [pos1 - pos2; pos3 - pos2], {sensor1; sensor3}); sensor3.plot(pos3(3), [pos1 - pos3; pos2 - pos3], {sensor1; sensor2}); end - function test_thesis_antenna_patterns(tc) - f = figure; - xi = -90:0.1:90; % degrees - n = [1, 4, 16, 64, 256]; - for ii = 1:length(n) - gain = 10 .* log10(2 * (n(ii) + 1)) + 10 .* n(ii) .* log10(max(0, cosd(xi))); - plot(xi, gain, "LineWidth", 2); - if ii == 1 - hold("on"); - end - end - hold("off"); - grid("on"); - ylim([-100, 35]); - xlim([-90, 90]); - xticks(-90:45:90); - % exportThesisFig switches all text to the LaTeX interpreter - title("Transmitting Antenna Gain Model"); - ylabel("Transmitting Antenna Gain (dBi)"); - xlabel("Elevation Angle from Antenna Boresight $\xi$ ($^\circ$)"); - legend("$\delta = " + string(n) + "$"); - - % exportThesisFig(f, "antennapatterns", 3/4); - % close(f); - end end end \ No newline at end of file diff --git a/test/test_sigmoidSensor.m b/test/test_sigmoidSensor.m index 6aeec87..523be69 100644 --- a/test/test_sigmoidSensor.m +++ b/test/test_sigmoidSensor.m @@ -57,30 +57,6 @@ classdef test_sigmoidSensor < matlab.unittest.TestCase % Performance at great tilt should be 0 tc.verifyEqual(tc.testClass.sensorPerformance([0, 0, h], [5, 5, 0]), 0, "AbsTol", 1e-9); end - function test_thesis_graphic(tc) - s = [sigmoidSensor, sigmoidSensor]; - - alphaDist = [2, 8]; - betaDist = [25, 1]; - alphaTilt = [50, 15]; % degrees - betaTilt = [0.2, 3]; - - s(1) = s(1).initialize(alphaDist(1), betaDist(1), alphaTilt(1), betaTilt(1)); - s(2) = s(2).initialize(alphaDist(2), betaDist(2), alphaTilt(2), betaTilt(2)); - - set(groot, 'defaultTextInterpreter', 'latex'); - set(groot, 'defaultAxesTickLabelInterpreter', 'latex'); - set(groot, 'defaultLegendInterpreter', 'latex'); - set(groot, 'defaultColorbarTickLabelInterpreter', 'latex'); - set(groot, 'defaultAxesFontSize', 10); % = \footnotesize in a 12pt doc - set(groot, 'defaultTextFontSize', 10); - - % Plot - f = s(1).plotParameters; - f = s(2).plotParameters(f); - - % f = exportThesisFig(f, "sigmoidmodel", 3/4); - end end end \ No newline at end of file diff --git a/test_thesis_plots.m b/test_thesis_plots.m new file mode 100644 index 0000000..4aa07a1 --- /dev/null +++ b/test_thesis_plots.m @@ -0,0 +1,230 @@ +classdef test_thesis_plots < matlab.unittest.TestCase + % test class that generates plots I use in my thesis + methods (Test) + function service_performance(tc) + % TX parameters + P_TX = 10e-3; + T_0 = 300; + BW = 20e6; + f_c = 2e9; + G_RX_dBi = 3; + beamwidthExponent = [16, 12, 32]; + elevation = [30, 20, 45]; + azimuth = [15, 325, 200]; + lossExponent = 2; + + % TX Positions + pos1 = [50, 50, 35]; + pos2 = [70, 30, 30]; + pos3 = [40, 75, 40]; + + % User Distribution Parameters + mu = [78, 66]; + sig = reshape([2000, 1100; 1100, 2000], [1, 2, 2]); + + sensor1 = rfSensor; + sensor1 = sensor1.initialize(P_TX, T_0, BW, f_c, G_RX_dBi, beamwidthExponent(1), elevation(1), azimuth(1), lossExponent); + sensor2 = rfSensor; + sensor2 = sensor2.initialize(P_TX, T_0, BW, f_c, G_RX_dBi, beamwidthExponent(2), elevation(2), azimuth(2), lossExponent); + sensor3 = rfSensor; + sensor3 = sensor3.initialize(P_TX, T_0, BW, f_c, G_RX_dBi, beamwidthExponent(3), elevation(3), azimuth(3), lossExponent); + + % Compute sensor and system performance for a single user distribution + domain = rectangularPrism; + domain = domain.initialize([zeros(1, 3); 100 * ones(1, 3)], "DOMAIN", "Domain"); + + objectiveFunction = objectiveFunctionWrapper(mu, sig); + domain.objective = domain.objective.initialize(objectiveFunction, domain, 0.1, 0, 1e-6, mu, sig); + + % Build target grid matching the domain objective grid and flatten user distribution + nX = size(domain.objective.values, 2); + nY = size(domain.objective.values, 1); + x_domain = linspace(0, 100, nX); + y_domain = linspace(0, 100, nY); + [Xd, Yd] = meshgrid(x_domain, y_domain); + targetPos_domain = [Xd(:), Yd(:), zeros(numel(Xd), 1)]; + userDist = domain.objective.values(:) / sum(domain.objective.values, "all"); + + [SINR1, ~, sensor1, others] = sensor1.sensorPerformance(pos1, targetPos_domain, [pos2; pos3], {sensor2; sensor3}); + sensor2 = others{1}; sensor3 = others{2}; + [SINR2, ~, sensor2, others] = sensor2.sensorPerformance(pos2, targetPos_domain, [pos1; pos3], {sensor1; sensor3}); + sensor1 = others{1}; sensor3 = others{2}; + [SINR3, ~, sensor3, ~] = sensor3.sensorPerformance(pos3, targetPos_domain, [pos1; pos2], {sensor1; sensor2}); + + maxSINR_linear = max([10.^(SINR1/10), 10.^(SINR2/10), 10.^(SINR3/10)], [], 2); + + % Plot SINR from each UAV's perspective with a shared color scale, + % one top-down map per UAV, in the style of the sigmoid sensor figure + sharedClim = [-80, max([SINR1; SINR2; SINR3])]; + SINR = {SINR1, SINR2, SINR3}; + sensors = {sensor1, sensor2, sensor3}; + positions = [pos1; pos2; pos3]; + for ii = 1:3 + f = plotSinrMap(x_domain, y_domain, reshape(SINR{ii}, size(Xd)), ... + sensors, positions, ii, sharedClim); + if ii == 3 % only plot for #3 right now + f = exportThesisFig(f, "sinrexample" + ii, 0.95, Colormap="parula", Recolor=false, Resolution=300); + end + close(f); + end + end + function sigmoid_sensor_performance(tc) + s = sigmoidSensor; + s = s.initialize(65, 1, 50, 0.2, 30, 120); + p = 50 * ones(1, 3); + mu = [50, 50]; + sig = 25 * reshape(eye(2), [1, 2, 2]); + domain = rectangularPrism; + domain = domain.initialize([zeros(1, 3); 100 * ones(1, 3)], "DOMAIN", "Domain"); + objectiveFunction = objectiveFunctionWrapper(mu, sig); + domain.objective = domain.objective.initialize(objectiveFunction, domain, 0.1, 0, 1e-6, mu, sig); + + S = s.sensorPerformance(p, [domain.objective.X(:), domain.objective.Y(:), zeros(size(domain.objective.X(:)))]); + S = reshape(S, size(domain.objective.X)); + + % meshgrid puts y down the rows, so plot against the grid vectors + % and put y = 0 at the bottom (imagesc defaults to row 1 at the top) + x = domain.objective.X(1, :); + y = domain.objective.Y(:, 1); + + f = figure; + imagesc(x, y, S); + set(gca, "YDir", "normal"); + axis("image"); + hold("on"); + % Boresight direction: the arrow runs from the UAV to where the + % boresight meets the ground, so its length grows with tilt + % (altitude * tan(tilt)). Azimuth 0 = +Y, 90 = +X, as in sensorPerformance. + reach = p(3) * tand(s.tilt); + hit = p(1:2) + reach * [sind(s.azimuth), cosd(s.azimuth)]; + quiver(p(1), p(2), hit(1) - p(1), hit(2) - p(2), 0, "k", "MaxHeadSize", 0.4); + hUav = plot(p(1), p(2), "k+"); + hHit = plot(hit(1), hit(2), "kx"); + hold("off"); + legend([hUav, hHit], "UAV position", "Boresight-ground intersection", "Location", "northwest"); + xlabel("$x$ (m)"); + ylabel("$y$ (m)"); + cb = colorbar; + cb.Label.String = "Sensor performance $s_n$"; + clim([0, 1]); + + % Full text width so it's included at 1:1 like the other figures; + % 0.8 fits the square axes plus colorbar without side margins + f = exportThesisFig(f, "sigmoidsensorexample", 0.8, Colormap="parula", Recolor=false, Resolution=300); + close(f); + end + function antenna_patterns(tc) + f = figure; + xi = -90:0.1:90; % degrees + n = [1, 4, 16, 64, 256]; + for ii = 1:length(n) + gain = 10 .* log10(2 * (n(ii) + 1)) + 10 .* n(ii) .* log10(max(0, cosd(xi))); + plot(xi, gain, "LineWidth", 2); + if ii == 1 + hold("on"); + end + end + hold("off"); + grid("on"); + ylim([-100, 35]); + xlim([-90, 90]); + xticks(-90:45:90); + % exportThesisFig switches all text to the LaTeX interpreter + title("Transmitting Antenna Gain Model"); + ylabel("Transmitting Antenna Gain (dBi)"); + xlabel("Elevation Angle from Antenna Boresight $\xi$ ($^\circ$)"); + legend("$\delta = " + string(n) + "$"); + + exportThesisFig(f, "antennapatterns", 3/4); + close(f); + end + function sigmoid_sensor_membership_functions(tc) + s = [sigmoidSensor, sigmoidSensor]; + + alphaDist = [2, 8]; + betaDist = [25, 1]; + alphaTilt = [50, 15]; % degrees + betaTilt = [0.2, 3]; + + s(1) = s(1).initialize(alphaDist(1), betaDist(1), alphaTilt(1), betaTilt(1)); + s(2) = s(2).initialize(alphaDist(2), betaDist(2), alphaTilt(2), betaTilt(2)); + + set(groot, 'defaultTextInterpreter', 'latex'); + set(groot, 'defaultAxesTickLabelInterpreter', 'latex'); + set(groot, 'defaultLegendInterpreter', 'latex'); + set(groot, 'defaultColorbarTickLabelInterpreter', 'latex'); + set(groot, 'defaultAxesFontSize', 10); % = \footnotesize in a 12pt doc + set(groot, 'defaultTextFontSize', 10); + + % Plot + f = s(1).plotParameters; + f = s(2).plotParameters(f); + + f = exportThesisFig(f, "sigmoidmodel", 3/4); + close(f); + end + end + +end +function f = plotSinrMap(x, y, SINR, sensors, positions, tx, cLimits) + % Top-down SINR map for transmitting UAV tx, matching the sigmoid sensor + % figure. Every UAV gets a marker, an arrow to the X where its boresight + % meets the ground, and its -3 dB (half-power) beam footprint; the + % transmitting UAV is drawn solid, the interferers dashed. + f = figure; + imagesc(x, y, SINR); + set(gca, "YDir", "normal"); % y = 0 at the bottom + axis("image"); + hold("on"); + for ii = 1:numel(sensors) + [hit, ring] = beamFootprint(sensors{ii}, positions(ii, :)); + pos = positions(ii, :); + if ii == tx + style = "-"; marker = "k+"; + else + style = "--"; marker = "ko"; + end + hRing(ii) = plot(ring(:, 1), ring(:, 2), "k" + style); %#ok + quiver(pos(1), pos(2), hit(1) - pos(1), hit(2) - pos(2), 0, "k", ... + "LineStyle", style, "MaxHeadSize", 0.4); + hHit = plot(hit(1), hit(2), "kx"); + hUav(ii) = plot(pos(1), pos(2), marker, "MarkerFaceColor", "w"); %#ok + end + hold("off"); + rx = setdiff(1:numel(sensors), tx); + legend([hUav(tx), hUav(rx(1)), hHit, hRing(tx), hRing(rx(1))], ... + "Transmitting UAV", "Interfering UAVs", "Boresight-ground intersection", ... + "Transmitter $-3$ dB contour", "Interferer $-3$ dB contour", ... + "Location", "southoutside", "NumColumns", 2); + xlim([min(x), max(x)]); % footprints can extend past the area + ylim([min(y), max(y)]); + xlabel("$x$ (m)"); + ylabel("$y$ (m)"); + cb = colorbar; + cb.Label.String = "SINR (dB)"; + clim(cLimits); +end + +function [hit, ring] = beamFootprint(sensor, pos) + % Where the boresight meets the ground (z = 0), and the ground curve where + % the gain is 3 dB below boresight: the intersection of the half-power cone + % with the ground, as in rfSensor/plot. Azimuth 0 = +Y, 90 = +X. + tlt = sensor.tilt; + az = sensor.azimuth; + hit = pos(1:2) + pos(3) * tand(tlt) * [sind(az), cosd(az)]; + + % Rotate nadir to the boresight + Ry = [cosd(tlt), 0, -sind(tlt); 0, 1, 0; sind(tlt), 0, cosd(tlt)]; + Rz = [sind(az), -cosd(az), 0; cosd(az), sind(az), 0; 0, 0, 1]; + R = Rz * Ry; + b = R * [0; 0; -1]; + u = R * [1; 0; 0]; + v = R * [0; 1; 0]; + + ha = sensor.halfAngle(); + phi = linspace(0, 2*pi, 720)'; + dirs = cosd(ha) .* b' + sind(ha) .* (cos(phi) .* u' + sin(phi) .* v'); + t = -pos(3) ./ dirs(:, 3); + t(t <= 0) = NaN; % cone generators that never reach the ground + ring = pos(1:2) + t .* dirs(:, 1:2); +end diff --git a/util/exportThesisFig.m b/util/exportThesisFig.m index 4e1dafe..970bbba 100644 --- a/util/exportThesisFig.m +++ b/util/exportThesisFig.m @@ -5,8 +5,12 @@ function f = exportThesisFig(f, name, aspect, opts) % exportThesisFig(f, "sigmoidmodel") % full width, 4:3 % exportThesisFig(f, "traj", 1, Width=230) % half-width subfigure % exportThesisFig(f, "traj", 3/4, LineStyles=true) % distinguishable in grayscale + % exportThesisFig(f, "perf", 0.8, Width=330) % imagesc plot -> name.png % - % The figure is modified in place. + % Figures containing images (imagesc, image) are exported as a lossless + % 400 dpi PNG instead: MATLAB's vector PDF export draws them with a + % visible diagonal seam. The PNG records its dpi, so LaTeX still places it + % at Width points. The figure is modified in place. arguments (Input) f (1, 1) {mustBeA(f, "matlab.ui.Figure")}; name (1, 1) string; @@ -15,11 +19,17 @@ function f = exportThesisFig(f, name, aspect, opts) opts.FontSize (1, 1) double {mustBePositive} = 10; % \footnotesize in a 12pt doc opts.LineWidth (1, 1) double {mustBePositive} = 2; % data lines, points opts.AxesLineWidth (1, 1) double {mustBePositive} = 0.75; - opts.MarkerSize (1, 1) double {mustBePositive} = 4; + opts.MarkerSize (1, 1) double {mustBePositive} = 4; % markers on lines + opts.PointMarkerSize (1, 1) double {mustBePositive} = 16; % marker-only series, e.g. plot(x, y, "k+") + opts.PointMarkerLineWidth (1, 1) double {mustBePositive} = 2.5; opts.Palette (:, 3) double = okabeIto(); % colorblind- and print-safe opts.Recolor (1, 1) logical = true; % also override explicitly set colors opts.LineStyles (1, 1) logical = false; % cycle line styles per series - opts.Raster (1, 1) logical = false; % for very dense plots + opts.Raster (1, 1) logical = ~isempty(findall(f, 'Type', 'image')); % PNG; also for very dense plots + % dpi, for Raster. MATLAB silently writes an all-black PNG when the image is + % too large (3909 px wide failed, 2750 px worked), so full width stays at 400 + opts.Resolution (1, 1) double {mustBePositive} = 400; + opts.Colormap = []; % for images/surfaces, e.g. "parula"; [] keeps the figure's opts.OutDir (1, 1) string = matlab.project.rootProject().RootFolder; end arguments (Output) @@ -41,6 +51,17 @@ function f = exportThesisFig(f, name, aspect, opts) if opts.LineStyles ax.LineStyleOrder = styles; end + if ~isempty(opts.Colormap) + colormap(ax, opts.Colormap); + end + end + + % Colorbars (image and surface plots) get the same frame as the axes. Their + % Label isn't reached by findall above, so set it directly. + for cb = findall(f, 'Type', 'colorbar')' + cb.LineWidth = opts.AxesLineWidth; + cb.Label.Interpreter = 'latex'; + cb.Label.FontSize = opts.FontSize; end % Data series. SeriesIndex is the series' position in its axes' color @@ -57,11 +78,21 @@ function f = exportThesisFig(f, name, aspect, opts) isLine = isprop(s, 'Color'); isScatter = isa(s, 'matlab.graphics.chart.primitive.Scatter'); - if isLine || isScatter - s.LineWidth = opts.LineWidth; - end - if isprop(s, 'MarkerSize') - s.MarkerSize = opts.MarkerSize; + % Marker-only series (no connecting line) mark individual points; a + % thick data-line width would fill in markers like + and x + isPoint = isLine && isprop(s, 'Marker') && strcmp(s.LineStyle, 'none') ... + && ~strcmp(s.Marker, 'none'); + + if isPoint + s.LineWidth = opts.PointMarkerLineWidth; + s.MarkerSize = opts.PointMarkerSize; + else + if isLine || isScatter + s.LineWidth = opts.LineWidth; + end + if isprop(s, 'MarkerSize') + s.MarkerSize = opts.MarkerSize; + end end if opts.LineStyles && isLine && ~strcmp(s.LineStyle, 'none') s.LineStyle = styles(mod(k - 1, numel(styles)) + 1); @@ -89,9 +120,9 @@ function f = exportThesisFig(f, name, aspect, opts) f.PaperSize = [w, h]; if opts.Raster - renderer = {'-image', '-r300'}; + % PNG rather than '-dpdf -image', which embeds a lossy JPEG + print(f, fullfile(opts.OutDir, name + ".png"), '-dpng', "-r" + opts.Resolution); else - renderer = {'-vector'}; + print(f, fullfile(opts.OutDir, name + ".pdf"), '-dpdf', '-vector'); end - print(f, fullfile(opts.OutDir, name + ".pdf"), '-dpdf', renderer{:}); end