Сжать набор точек многоугольника в более короткий набор точек

У двух существующих ответов есть большие недостатки:

• Метод Дурки работает только в том случае, если расстояния между последующими точками точно идентичны. Точки должны иметь координаты, идеально представляемые в виде значений с плавающей запятой, чтобы расстояния между последующими точками на линии могли быть идентичными. Если точки не равноудалены, метод ничего не делает. Кроме того, начало и конец многоугольника не исследуются вместе, поэтому, если прямая линия образуется через начало / конец многоугольника, останется на одну точку слишком много.

• Метод ShadowMan лучше тем, что расстояния не обязательно должны быть идентичными, а линия, пересекающая начало / конец многоугольника, является правильно обработан. Однако он также использует сравнения на равенство с плавающей запятой, что в общем случае не работает. Только с целочисленными координатами этот метод будет работать правильно. Кроме того, он использует vecnorm (который вычисляет квадратный корень) и деление, оба являются относительно дорогостоящими операциями (по сравнению с методом, показанным здесь).

Чтобы увидеть, образуют ли три точки прямую линию, можно использовать простое арифметическое правило. Допустим, у нас есть точки p0, p1 и p2. Вектор от p0 до p1 и вектор от p0 до p2 образуют основу параллелограмма, площадь которого можно вычислить с помощью перекрестное произведение двух векторов (в 2D подразумевается, что в перекрестном произведении используется z=0, с результирующим вектором, имеющим x=0 и y=0, полезно только значение z; таким образом, 2D Предполагается, что перекрестное произведение дает скалярное значение). Его можно вычислить следующим образом:

v1 = p1 - p0;
v2 = p2 - p0;
x = v1(1)*v2(2) - v1(2)*v2(1);

x, перекрестное произведение, будет равно нулю, если два вектора параллельны, что означает, что три точки лежат на одной прямой. Но проверка на равенство нулю должна иметь некоторый допуск, поскольку арифметика с плавающей запятой неточна. Здесь я использую 1e-6 как допуск. Используйте значение, которое на несколько порядков меньше расстояния между вашими точками.

Учитывая входной набор точек p, мы можем найти угловые точки с помощью:

p1 = p;                                  % point 1
p0 = circshift(p1,1);                    % point 0
v1 = p1 - p0;                            % vector from point 0 to 1
v2 = circshift(p1,-1) - p0;              % vector from point 0 to 2
x = v1(:,1).*v2(:,2) - v1(:,2).*v2(:,1); % cross product
idx = abs(x) > 1e-6;                     % comparison with tolerance
p = p(idx,:);                            % corner points

Обратите внимание, что этот тест перекрестного произведения не удастся, если две последовательные точки имеют одинаковые координаты (т.е. один из векторов имеет нулевую длину). Если в данных могут быть дублированные точки, потребуется дополнительный тест.

Вот результаты трех методов. Я создал многоугольник с нетривиальными координатами и неравномерно расположенными вершинами. Я также поставил зазор между началом и концом в середине прямой кромки. Эти характеристики призваны показать недостатки двух других методов.

Это код, который я использовал для построения графика:

% Make a polygon that will be difficult for the other two methods
p = [0,0 ; 0.5,0 ; 1,0 ; 1,1 ; 0.5,1 ; 0,1];
p = p + rand(size(p))/3;
p(end+1,:) = p(1,:);
q = [];
for ii = 1:size(p,1)-1
   t = p(ii,:) + (p(ii+1,:) - p(ii,:)) .* [0;0.1;1/3;0.45;0.5897545;pi/4;exp(1)/3];
   q = [q;t];
end
q = circshift(q,3,1);

figure
subplot(2,2,1)
plot(q(:,1),q(:,2),'bo-')
axis equal
title('input')

subplot(2,2,2)
res1 = method1(q);
plot(res1(:,1),res1(:,2),'ro-')
axis equal
title('Durkee''s method')

subplot(2,2,3)
res2 = method2(q);
plot(res2(:,1),res2(:,2),'ro-')
axis equal
title('ShadowMan''s method')

subplot(2,2,4)
res3 = method3(q);
plot(res3(:,1),res3(:,2),'go-')
axis equal
title('correct method')

% Durkee's method: https://stackoverflow.com/a/55603145/7328782
function P = method1(P)
a = logical([1 diff(P(:,1),2)' 1]);
b = logical([1 diff(P(:,2),2)' 1]);
idx = or(a,b);
P = P(idx,:);
end

% ShadowMan's method: https://stackoverflow.com/a/55603040/7328782
function corners = method2(poly0)
poly0Z = circshift(poly0,1);
poly0I = circshift(poly0,-1);
unitVectIn =(poly0 - poly0I)./vecnorm((poly0 - poly0I),2,2);
unitVectOut =(poly0Z - poly0)./vecnorm((poly0Z - poly0),2,2);
cornerIndices = sum(unitVectIn == unitVectOut,2)==0;
corners = poly0(cornerIndices,:);
end
% vecnorm is new to R2017b, I'm still running R2017a.
function p = vecnorm(p,n,d)
% n is always 2
p = sqrt(sum(p.^2,d));
end

function p = method3(p1)
p0 = circshift(p1,1);
v1 = p1 - p0;
v2 = circshift(p1,-1) - p0;
x = v1(:,1).*v2(:,2) - v1(:,2).*v2(:,1);
idx = abs(x) > 1e-6;
p = p1(idx,:);
end

Хорошо, я адаптировал это для работы с неквадратными углами.

Рассмотрим треугольник, обозначенный точками

P = [0 0; 1 0; 2 0; 1.5 1; 1 2; .5 1; 0 0];

Это массив 7x2, если мы затем определим 2 производных вектора, как определено в вопросе, который я упомянул в комментариях.

a = logical([1 diff(P(:,1),2)' 1]);
b = logical([1 diff(P(:,2),2)' 1]);

Оттуда мы можем объединить их, чтобы получить новую переменную индексации

idx = or(a,b);

Наконец, мы можем использовать это, чтобы сделать наш сюжет

line(P(idx,1), P(idx,2),'Color','k','LineWidth',3,'LineStyle',':');

Если вы делаете линейный график, я думаю, вам нужно установить для последней переменной значение false.

idx(end) = false;

«Векторный» способ можно сделать довольно изящно. Я тоже пробовал использовать цикл for, и вы можете сделать то же самое, но вы попросили вектор, так что вот мой способ сделать это.

Единственное изменение, которое я внес в ваши данные, - это удалить все реплики перед запуском этого скрипта. Кроме того, указанные точки должны быть расположены в последовательном порядке по часовой стрелке или против часовой стрелки.

    poly0Z = circshift(poly0,1);
    poly0I = circshift(poly0,-1); 
    unitVectIn =(poly0 - poly0I)./vecnorm((poly0 - poly0I),2,2);
    unitVectOut =(poly0Z - poly0)./vecnorm((poly0Z - poly0),2,2)  ;
    cornerIndices = sum(unitVectIn == unitVectOut,2)==0
    corners = poly0(cornerIndices,:)

    line(poly0(:,1), poly0(:,2),'Color','k','LineWidth',2,'LineStyle',':');
    hold on
    scatter(corners(:,1), corners(:,2),'filled')

В основе этого метода лежит переход к каждой точке, вычисление входящего единичного вектора и выходящего единичного вектора. Точки, в которых единичный вектор in не совпадает с единичным вектором out, являются углами.