cache RSS data for efficiency in computing all timestep SINRs

@rfSensor/rfSensor.m

@@ -14,11 +14,13 @@ classdef rfSensor
         % Values computed at initialization
         P_TX_dBm = NaN; % Transmit power (dBm)
         N = NaN; % Thermal noise
+        % Cached state (per timestep)
+        rssCache (:,1) double = double.empty(0,1); % linear-scale RSS to last target grid
     end
 
     methods (Access = public)
         [obj] = initialize(obj, txPower, bandwidth, centerFreq, rxGain); % initialize sensor, define parameters
-        [SINR, SNR] = sensorPerformance(obj, agentPos, targetPos, otherSensorsPos, otherSensors); % determine sensor performance for a given single sensor and target geometry
+        [SINR, SNR, obj, otherSensors] = sensorPerformance(obj, agentPos, targetPos, otherSensorsPos, otherSensors); % determine sensor performance for a given single sensor and target geometry
         [d, t, a] = computePointToPoints(obj, agentPos, targetPos);
         [f] = plotParameters(obj); % debug, plot sensor response as a function of distance and tilt angle
         [f] = plotPerformance(obj, altitude, otherSensorsPos, otherSensors); % debug, plot SNR or SINR ground heatmap for a given geometry

@rfSensor/sensorPerformance.m

@@ -1,4 +1,4 @@
-function [SINR, SNR] = sensorPerformance(obj, agentPos, targetPos, otherSensorsPos, otherSensors)
+function [SINR, SNR, obj, otherSensors] = sensorPerformance(obj, agentPos, targetPos, otherSensorsPos, otherSensors)
     arguments (Input)
         obj (1, 1) {mustBeA(obj, "rfSensor")};
         agentPos (1, 3) double;
@@ -9,19 +9,25 @@ function [SINR, SNR] = sensorPerformance(obj, agentPos, targetPos, otherSensorsP
     arguments (Output)
         SINR (:, 1) double;
         SNR (:, 1) double;
+        obj (1, 1) {mustBeA(obj, "rfSensor")};
+        otherSensors (:, 1) cell;
     end
     assert(size(otherSensorsPos, 1) == size(otherSensors, 1), "Mismatch in number of other sensor positions (%d) and number of other sensors (%d) provided", size(otherSensorsPos, 1), size(otherSensors, 1));
 
     [d, t, a] = obj.computePointToPoints(agentPos, targetPos);
 
-    % Performance is measured as SINR for this sensor
+    if isempty(obj.rssCache)
-    %% TODO: how should interference calculation be modified for
+        obj.rssCache = 10 .^ (0.1 .* obj.RSS(d, t, a));
-    % interference sources with different center frequencies and bandwidths?
+    end
-    S = 10 .^ (0.1 .* obj.RSS(d, t, a)); % Signal
+    S = obj.rssCache;
-    I = zeros(size(d)); % Interference from other agents
+
+    I = zeros(size(d));
     for ii = 1:size(otherSensors, 1)
-        [d_other, t_other, a_other] = otherSensors{ii}.computePointToPoints(otherSensorsPos(ii, 1:3), targetPos);
+        if isempty(otherSensors{ii}.rssCache)
-        I = I + 10 .^ (0.1 .* otherSensors{ii}.RSS(d_other, t_other, a_other));
+            [d_other, t_other, a_other] = otherSensors{ii}.computePointToPoints(otherSensorsPos(ii, 1:3), targetPos);
+            otherSensors{ii}.rssCache = 10 .^ (0.1 .* otherSensors{ii}.RSS(d_other, t_other, a_other));
+        end
+        I = I + otherSensors{ii}.rssCache;
     end
 
     SINR = 10*log10(S ./ (I + obj.N));

test/test_rfSensor.m

@@ -74,5 +74,74 @@ classdef test_rfSensor < matlab.unittest.TestCase
 
             tc.testClass.plotPerformance(altitude, otherSensorsPos, otherSensors);
         end
+
+        function plot_SINR_heterogenous_interferers_efficiently(tc)
+            P_TX = 1e-3;
+            BW = 20e6;
+            f_c = 2e9;
+            G_RX_dBi = 3;
+            altitude = 30;
+
+            sensor1 = rfSensor;
+            sensor1 = sensor1.initialize(P_TX, BW, f_c, G_RX_dBi, 0, 0);
+            sensor2 = rfSensor;
+            sensor2 = sensor2.initialize(P_TX, BW, f_c, G_RX_dBi, 0, 0);
+            sensor3 = rfSensor;
+            sensor3 = sensor3.initialize(P_TX, BW, f_c, G_RX_dBi, 0, 0);
+
+            pos1 = [0,  0,  altitude];
+            pos2 = [6, -4,  altitude - 1];
+            pos3 = [-2, 6,  altitude];
+
+            % Build a shared target grid
+            distances = -15:0.25:15;
+            [Xg, Yg] = meshgrid(distances, distances);
+            targetPos = [Xg(:), Yg(:), zeros(numel(Xg), 1)];
+
+            % Call 1: cache empty, does all computations for this timestep
+            [SINR1, ~, sensor1, others] = sensor1.sensorPerformance(pos1, targetPos, [pos2; pos3], {sensor2; sensor3});
+            sensor2 = others{1};
+            sensor3 = others{2};
+
+            % Calls 2 and 3 use cached data
+            [SINR2, ~, sensor2, others] = sensor2.sensorPerformance(pos2, targetPos, [pos1; pos3], {sensor1; sensor3});
+            sensor1 = others{1};
+            sensor3 = others{2};
+
+            [SINR3, ~, sensor3, ~] = sensor3.sensorPerformance(pos3, targetPos, [pos1; pos2], {sensor1; sensor2});
+
+
+            % All caches should be populated after the three calls
+            tc.assertNotEmpty(sensor1.rssCache);
+            tc.assertNotEmpty(sensor2.rssCache);
+            tc.assertNotEmpty(sensor3.rssCache);
+
+            % Plot SINR from each UAV's perspective
+            sz = size(Xg);
+            SINR1 = reshape(SINR1, sz);
+            SINR2 = reshape(SINR2, sz);
+            SINR3 = reshape(SINR3, sz);
+
+            f = figure;
+            tiledlayout(f, 1, 3, TileSpacing="compact", Padding="compact");
+
+            nexttile;
+            imagesc(distances, distances, SINR1); axis image; set(gca, YDir="normal"); hold on;
+            scatter(pos1(1), pos1(2), 80, "g", "o", LineWidth=2);
+            scatter([pos2(1), pos3(1)], [pos2(2), pos3(2)], 80, "r", "x", LineWidth=2);
+            hold off; cb = colorbar; cb.Label.String = "SINR (dB)"; xlabel("X (m)"); ylabel("Y (m)"); title("SINR: UAV 1");
+
+            nexttile;
+            imagesc(distances, distances, SINR2); axis image; set(gca, YDir="normal"); hold on;
+            scatter(pos2(1), pos2(2), 80, "g", "o", LineWidth=2);
+            scatter([pos1(1), pos3(1)], [pos1(2), pos3(2)], 80, "r", "x", LineWidth=2);
+            hold off; cb = colorbar; cb.Label.String = "SINR (dB)"; xlabel("X (m)"); ylabel("Y (m)"); title("SINR: UAV 2");
+
+            nexttile;
+            imagesc(distances, distances, SINR3); axis image; set(gca, YDir="normal"); hold on;
+            scatter(pos3(1), pos3(2), 80, "g", "o", LineWidth=2);
+            scatter([pos1(1), pos2(1)], [pos1(2), pos2(2)], 80, "r", "x", LineWidth=2);
+            hold off; cb = colorbar; cb.Label.String = "SINR (dB)"; xlabel("X (m)"); ylabel("Y (m)"); title("SINR: UAV 3");
+        end
     end
 end