Indexing the intersection of polygons.

I have the following code that creates two polygons and then intersection. I want to do indexing such that I have the index of intersection points in any of the input polygons.

%%
clc
clear all
close all
%%
xlimit = [3 13];
ylimit = [2  8];
xbox = xlimit([1 1 2 2 1]);
ybox = ylimit([1 2 2 1 1]);

distancesA = [0; cumsum(sqrt(diff(xbox').^2 + diff(ybox').^2))];
numPointsA = distancesA(end)*100*2;
equalSpacedDistancesA = linspace(0, distancesA(end), numPointsA);
xEqualSpacedA = interp1(distancesA, xbox, equalSpacedDistancesA);
yEqualSpacedA = interp1(distancesA, ybox, equalSpacedDistancesA);

hold on
plot(xEqualSpacedA, yEqualSpacedA, 'DisplayName','poly1','Marker','.')

x = [0 6  4  8 8 10 14 10 14 4 4 6 9 15];
y = [4 6 10 11 7  6 10 10  6 0 3 4 3  6];

distancesB = [0; cumsum(sqrt(diff(x').^2 + diff(y').^2))];
numPointsB = distancesB(end)*100*2;
equalSpacedDistancesB = linspace(0, distancesB(end), numPointsB);
xEqualSpacedB = interp1(distancesB, x, equalSpacedDistancesB);
yEqualSpacedB = interp1(distancesB, y, equalSpacedDistancesB);

hold on
plot(xEqualSpacedB, yEqualSpacedB, 'DisplayName','poly2','Marker','.')

%
[xi,yi] = polyxpoly(xEqualSpacedA,yEqualSpacedA,xEqualSpacedB,yEqualSpacedB);
scatter(xi,yi,'Marker','o')

Answer 1

[Hitesh]

Hi @Prasanna Routray,

You need to calculate the "Euclidean distance" from each intersection point to every point in the first polygon (xEqualSpacedA, yEqualSpacedA). Similarly, perform this calculation for the second polygon. After computing these distances, you need to use the "min" function to find the index of the nearest point in each polygon to each intersection point. This approach will help you to accurately identify the location of intersection points relative to the vertices of the polygons.Kindly refer to the revised code for finding indexing of intersection polygons.

clc;
clear all;
close all;

% Define the limits and create the first polygon
xlimit = [3 13];
ylimit = [2  8];
xbox = xlimit([1 1 2 2 1]);
ybox = ylimit([1 2 2 1 1]);
distancesA = [0; cumsum(sqrt(diff(xbox').^2 + diff(ybox').^2))];
numPointsA = distancesA(end)*100*2;
equalSpacedDistancesA = linspace(0, distancesA(end), numPointsA);
xEqualSpacedA = interp1(distancesA, xbox, equalSpacedDistancesA);
yEqualSpacedA = interp1(distancesA, ybox, equalSpacedDistancesA);

% Plot the first polygon
hold on;
plot(xEqualSpacedA, yEqualSpacedA, 'DisplayName','poly1','Marker','.');

% Define and create the second polygon
x = [0 6  4  8 8 10 14 10 14 4 4 6 9 15];
y = [4 6 10 11 7  6 10 10  6 0 3 4 3  6];
distancesB = [0; cumsum(sqrt(diff(x').^2 + diff(y').^2))];
numPointsB = distancesB(end)*100*2;
equalSpacedDistancesB = linspace(0, distancesB(end), numPointsB);
xEqualSpacedB = interp1(distancesB, x, equalSpacedDistancesB);
yEqualSpacedB = interp1(distancesB, y, equalSpacedDistancesB);

% Plot the second polygon
hold on;
plot(xEqualSpacedB, yEqualSpacedB, 'DisplayName','poly2','Marker','.');

% Calculate intersection points
[xi, yi] = polyxpoly(xEqualSpacedA, yEqualSpacedA, xEqualSpacedB, yEqualSpacedB);

% Scatter plot intersection points
scatter(xi, yi, 'Marker', 'o');

% Find indices of intersection points in the original polygons
indicesA = [];
indicesB = [];
for k = 1:length(xi)
    % Find the closest point in polygon A
    distancesToA = sqrt((xEqualSpacedA - xi(k)).^2 + (yEqualSpacedA - yi(k)).^2);
    [~, indexA] = min(distancesToA);
    indicesA = [indicesA; indexA];
    
    % Find the closest point in polygon B
    distancesToB = sqrt((xEqualSpacedB - xi(k)).^2 + (yEqualSpacedB - yi(k)).^2);
    [~, indexB] = min(distancesToB);
    indicesB = [indicesB; indexB];
end

% Display the indices
disp('Indices in Polygon A:');
disp(indicesA);
disp('Indices in Polygon B:');
disp(indicesB);

Answer 2

[KALYAN ACHARJYA]

Other way:

xlimit = [3 13]; ylimit = [2 8]; xbox = xlimit([1 1 2 2 1]); ybox = ylimit([1 2 2 1 1]);

distancesA = [0; cumsum(sqrt(diff(xbox').^2 + diff(ybox').^2))]; 
numPointsA = distancesA(end)*100*2;
equalSpacedDistancesA = linspace(0, distancesA(end), numPointsA);
xEqualSpacedA = interp1(distancesA, xbox, equalSpacedDistancesA);
yEqualSpacedA = interp1(distancesA, ybox, equalSpacedDistancesA);

plot(xEqualSpacedA, yEqualSpacedA, 'DisplayName','poly1','Marker','.')
hold on;
x = [0 6 4 8 8 10 14 10 14 4 4 6 9 15]; 
y = [4 6 10 11 7 6 10 10 6 0 3 4 3 6];

distancesB = [0; cumsum(sqrt(diff(x').^2 + diff(y').^2))]; 
numPointsB = distancesB(end)*100*2; 
equalSpacedDistancesB = linspace(0, distancesB(end), numPointsB);
xEqualSpacedB = interp1(distancesB, x, equalSpacedDistancesB);
yEqualSpacedB = interp1(distancesB, y, equalSpacedDistancesB);
plot(xEqualSpacedB, yEqualSpacedB, 'DisplayName','poly2','Marker','.')

% [xi,yi] = polyxpoly(xEqualSpacedA,yEqualSpacedA,xEqualSpacedB,yEqualSpacedB); scatter(xi,yi,'Marker','o')

%% 
% Find intersection points
[xi, yi] = polyxpoly(xEqualSpacedA, yEqualSpacedA, xEqualSpacedB, yEqualSpacedB);
distancesA = sqrt((xEqualSpacedA(:) - xi').^2 + (yEqualSpacedA(:) - yi').^2); % Ensure column vectors
[~, indicesA] = min(distancesA, [], 1);  % Find closest points
distancesB = sqrt((xEqualSpacedB(:) - xi').^2 + (yEqualSpacedB(:) - yi').^2); % Ensure column vectors
[~, indicesB] = min(distancesB, [], 1);  % Find closest points

% Make sure that the intersection points and their indices match in size
numIntersections = length(xi); 
indicesA = indicesA(:); 
indicesB = indicesB(:); 

% This happen depending on data: If there are fewer indices than intersections, pad with NaN
if length(indicesA) < numIntersections
    indicesA = [indicesA; NaN(numIntersections - length(indicesA), 1)];
end
if length(indicesB) < numIntersections
    indicesB = [indicesB; NaN(numIntersections - length(indicesB), 1)];
end

% Display the results (If it Required, otherwise ignore it)
disp('Indices of intersection points:');
disp(table(xi, yi, indicesA, indicesB, 'VariableNames', {'X_Intersect', 'Y_Intersect', 'Index_A', 'Index_B'}));

% Show intersection points
scatter(xi, yi, 50,'MarkerEdgeColor',[0 .5 .5],'MarkerFaceColor',[0 .7 .7],'LineWidth',1.5);