Find only interior points from set of points

14 vues (au cours des 30 derniers jours)
Ryan
Ryan le 18 Avr 2023
Commenté : Image Analyst le 18 Avr 2023
I have a collection of interior and exterior points and need to remove all exterior points from the set. The points can be in any order and any shape. The exterior shape does not have to match the interior shape. However, the inner shape NEVER intersects the outer shape. Both shapes share a common "center" point and I have that data. For example, I could have the points like...
t1 = linspace(0,2*pi, 100);
t2 = wrapTo2Pi(t1 + 0.05);
r1 = 5;
r2 = 10;
r3 = 7;
r4 = 12;
x1 = r1*cos(t1);
y1 = r2*sin(t1);
x2 = r3*cos(t2);
y2 = r4*sin(t2);
xGiven = [x1, x2];
yGiven = [y1, y2];
figure()
scatter(xGiven , yGiven);
Assume, that the only information provided is xGiven, and yGiven. No information on how they were constructed is provided.
Here, I only want the points of the interior ring. The common point shared is [0,0] since both shapes are symmetric about that point.
The only toolbox I have access to is the signal processing toolbox.
Any suggestions would be greatly appreciated.

Réponse acceptée

Image Analyst
Image Analyst le 18 Avr 2023
Try this:
% Demo by Image Analyst separate point sets using the convex hull.
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format short g;
format compact;
fontSize = 16;
markerSize = 16;
t1 = linspace(0,2*pi, 100);
t2 = wrapTo2Pi(t1 + 0.05);
r1 = 5;
r2 = 10;
r3 = 7;
r4 = 12;
x1 = r1*cos(t1);
y1 = r2*sin(t1);
x2 = r3*cos(t2);
y2 = r4*sin(t2);
xGiven = [x1, x2];
yGiven = [y1, y2];
figure()
subplot(3, 1, 1);
scatter(xGiven , yGiven, 'filled');
xlim([-8, 8]);
title('All Points', 'FontSize', fontSize);
grid on;
% First get exterior points. They are the convex hull.
chIndexes = convhull(xGiven , yGiven);
xExterior = xGiven(chIndexes);
yExterior = yGiven(chIndexes);
% Now get the interior points. They are whatever is not the convex hull
interiorIndexes = setdiff(1:length(xGiven), chIndexes);
xInterior = xGiven(interiorIndexes);
yInterior = yGiven(interiorIndexes);
% Plot red dots for the exterior points
subplot(3, 1, 2);
plot(xExterior, yExterior, 'r.', 'MarkerSize', markerSize);
xlim([-8, 8]);
grid on;
title('Exterior Points', 'FontSize', fontSize);
% Plot magenta dots for the interior points.
subplot(3, 1, 3);
plot(xInterior, yInterior, 'm.', 'MarkerSize', markerSize);
xlim([-8, 8]);
grid on;
title('Interior Points', 'FontSize', fontSize);
  2 commentaires
Ryan
Ryan le 18 Avr 2023
Thank you Image Analyst, is there a reason that there is one stray point in the interior points plot you produced?
Image Analyst
Image Analyst le 18 Avr 2023
I was wondering that myself, but I didn't have time to delve into it. I'd have to zoom way in to see if that point is actually not on the convex hull, or if it is, if there is a bug in the function.

Connectez-vous pour commenter.

Plus de réponses (1)

Akira Agata
Akira Agata le 18 Avr 2023
As an another possible solution, how about using boundary function to detect exterior points?
The following is an example:
% Sample data
rng('default');
x = rand(100, 1);
y = rand(100, 1);
% Identiry exterior points
k = boundary(x, y);
% Create index
idxExterior = false(size(x));
idxExterior(k) = true;
idxInterior = ~idxExterior;
% Show the result
figure
tiledlayout('flow')
ax1 = nexttile;
scatter(x,y, 'b.')
daspect([1 1 1])
title('All points')
ax2 = nexttile;
scatter(x(idxExterior), y(idxExterior), 'r.')
daspect([1 1 1])
title('Exterior points')
ax3 = nexttile;
scatter(x(idxInterior), y(idxInterior), 'm.')
daspect([1 1 1])
title('Interior points')
linkaxes([ax1, ax2, ax3], 'xy')
  2 commentaires
Ryan
Ryan le 18 Avr 2023
What is the difference between using convex hull and boundary?
Image Analyst
Image Analyst le 18 Avr 2023
Looks like with boundary the points don't have to be convex. The "outer" points can go in and out, have protrusions and bays.
help boundary
BOUNDARY Boundary of a set of points in 2D/3D space K = BOUNDARY(X,Y) returns a single conforming boundary around the points (X,Y). X and Y are column vectors of equal size that specify the coordinates. K is a vector of point indices that represents a compact boundary around the points. Unlike the convex hull, the boundary can shrink towards the interior of the hull to envelop the points. K = BOUNDARY(X,Y,Z) returns a single conforming boundary around the points (X,Y,Z). X, Y and Z are column vectors of equal size that specify the coordinates. K is a triangulation that represents a compact boundary around the points. K is of size mtri-by-3, where mtri is the number of triangular facets. That is, each row of K is a triangle defined in terms of the point indices. Unlike the convex hull, the boundary can shrink towards the interior of the hull to envelop the points. K = BOUNDARY(P) returns the 2D/3D boundary of the points P. This syntax is equivalent to the (X,Y) and (X,Y,Z) syntaxes where the columns of P are X,Y, or X,Y,Z respectively. K = BOUNDARY(...,S) provides an option of specifying the shrink factor S. The scalar S has a value in the range 0<=S<=1. Setting S to 0 gives the convex hull, while setting S to 1 gives a compact boundary that envelops the points. The default shrink factor is 0.5. [K,V] = BOUNDARY(...) returns the boundary K and the corresponding area/volume V bounded by K. Example 1: x = rand(30,1); y = rand(30,1); k = boundary(x,y); plot(x,y, '.', x(k), y(k), '-r') axis equal % Repeat this example by separately passing % a shrink factor of 0 and then 1 to the % boundary function. Example 2: x = rand(25,1); y = rand(25,1); z = rand(25,1); k = boundary(x,y,z,0.8); trisurf(k,x,y,z, 'Facecolor','cyan','FaceAlpha',0.8); axis equal; hold on plot3(x,y,z,'.r') See also alphaShape, triangulation, delaunayTriangulation, TRISURF, CONVHULL Documentation for boundary doc boundary Other uses of boundary paretotails/boundary polyshape/boundary piecewisedistribution/boundary

Connectez-vous pour commenter.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by