Looking for where line intersects on itself
显示 更早的评论
I am looking to see where a line that is plotted by x and y intersects on itself. I have tried using unique() with the index given to see where the same combination of x and y occur, but based on the data that index did not make sense. I have also tried linexline() which was for two lines intersecting not a line intersecting itself. I have also tried to just zoom in on the graph and find the point. Although the amount of detail required for that was not working on my graph.
4 个评论
Caitlin Bemis
2022-7-13
Star Strider
2022-7-15
This answer was posted while I was sleeping, and since I cannot improve on it, I deleted my answer.
Caitlin Bemis
2022-7-18
采纳的回答
If you interpolate your data by line segments it cross it selft at least 135 times
data = readtable('data.csv');
xy=table2array(data)';
X=polyselfx(xy);
fprintf('Number of intersections = %d\n', size(X,2))
figure
plot(xy(1,:),xy(2,:))
hold on
plot(X(1,:),X(2,:),'rx')
function X = polyselfx(P)
% Empty buffer
X = zeros(2,0);
filled = 0;
sizec = 0;
for n=2:size(P,2)-1
cn = seg2poly(P(:,n:n+1), P(:,1:n-1));
m = size(cn,2);
filled = filled+m;
% Buffer too small
if sizec < filled
sizec = max(filled, 2*sizec);
X(2,sizec) = 0;
end
% Store the result
X(:,filled+(-m+1:0)) = cn;
end
% remove the tails
X(:,filled+1:end) = [];
end % polyselfx
%%
function [cross_X, cross_t1, cross_t2, cross_i2] = seg2poly(s1, P)
a = s1(:,1);
M = P-a;
b = s1(:,2)-a;
nb2 = b(1)^2+b(2)^2;
% Check if the points are on the left/right side
x = [b(2) -b(1)]*M; % 1 x n
sx = sign(x);
% x -coordinates has opposite signs
crossx = sx(1:end-1).*sx(2:end) <= 0;
ix = find(crossx);
% cross point to the y-axis (along the segment)
x1 = x(ix);
x2 = x(ix+1);
d = b.'/nb2;
y1 = d*M(:,ix);
y2 = d*M(:,ix+1);
dx = x2-x1;
t1 = (y1.*x2-y2.*x1)./dx;
% Check if the cross point is inside the segment
ind = t1>0 & t1<=1;
if any(ind)
cross_t1 = t1(ind);
cross_t2= -x1(ind)./dx(ind);
cross_X = a + b*cross_t1;
cross_i2 = ix(ind);
else
cross_X = zeros(2,0);
cross_t1 = zeros(1,0);
cross_t2 = zeros(1,0);
cross_i2 = [];
end
end % seg2poly
5 个评论
Caitlin Bemis
2022-7-17
Caitlin Bemis
2022-7-17
Bruno Luong
2022-7-17
Inside the loop when cn is not impty the indexes are:
- n and cross_i2, the 4th output argument of seg2poly which is not retreived in my code. If you need it just put them in array.
- cross_t1 and cross_t2 are the fractional parts of the index where the intersection occurs.
Caitlin Bemis
2022-7-18
Bruno Luong
2022-7-18
编辑：Bruno Luong
2022-7-18
Here we go I modify the code and you get in idx the (first) indexes of the line segments and the fractional parametriec coordinates t of the intersections
data = readtable('data.csv');
xy=table2array(data)';
[X,idx,t]=polyselfx(xy);
fprintf('Number of intersections = %d\n', size(X,2))
figure
plot(xy(1,:),xy(2,:))
hold on
plot(X(1,:),X(2,:),'rx')
function [X, loc, t] = polyselfx(P)
% Crossing points polygons P
%
% INPUTS:
% P: two-row arrays, each column is a vertices
% OUTPUTS:
% X is two-row array, each column is an crossing point
% loc: two-row array, which edges the crossing point belong?
% first row corresponds to P1, second row to P2
% edge#1 is P(:,[1 2]), edge#2 is P(:,[2 3]), ... etc
% t: floating parametric of the crossing point
% Empty buffer
X = zeros(2,0);
loc = zeros(2,0);
t = zeros(2,0);
filled = 0;
sizec = 0;
for n=2:size(P,2)-1
[cn, t1, t2, i2] = seg2poly(P(:,n:n+1), P(:,1:n-1));
m = size(cn,2);
filled = filled+m;
% Buffer too small
if sizec < filled
sizec = max(filled, 2*sizec);
X(2,sizec) = 0;
loc(2,sizec) = 0;
t(2,sizec) = 0;
end
% Store the result
ifill = filled+(-m+1:0);
X(:,ifill) = cn;
loc(1,ifill) = n;
loc(2,ifill) = i2;
t(1,ifill) = t1;
t(2,ifill) = t2;
end
% remove the tails
X(:,filled+1:end) = [];
loc(:,filled+1:end) = [];
t(:,filled+1:end) = [];
end % polyselfx
%%
function [cross_X, cross_t1, cross_t2, cross_i2] = seg2poly(s1, P)
a = s1(:,1);
M = P-a;
b = s1(:,2)-a;
nb2 = b(1)^2+b(2)^2;
% Check if the points are on the left/right side
x = [b(2) -b(1)]*M; % 1 x n
sx = sign(x);
% x -coordinates has opposite signs
crossx = sx(1:end-1).*sx(2:end) <= 0;
ix = find(crossx);
% cross point to the y-axis (along the segment)
x1 = x(ix);
x2 = x(ix+1);
d = b.'/nb2;
y1 = d*M(:,ix);
y2 = d*M(:,ix+1);
dx = x2-x1;
t1 = (y1.*x2-y2.*x1)./dx;
% Check if the cross point is inside the segment
ind = t1>0 & t1<=1;
if any(ind)
cross_t1 = t1(ind);
cross_t2= -x1(ind)./dx(ind);
cross_X = bsxfun(@plus, a, b*cross_t1);
cross_i2 = ix(ind);
else
cross_X = zeros(2,0);
cross_t1 = zeros(1,0);
cross_t2 = zeros(1,0);
cross_i2 = [];
end
end % seg2poly
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Matrix Indexing 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
