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.

[Debadipto]

Please refer to this answer to a similar question.

[Caitlin Bemis]

This does answer the question partly, but I also need the index to create new subsets of data.

My end goal is to find the AUC of the new subset of data to calculate work.

[Star Strider]

This answer was posted while I was sleeping, and since I cannot improve on it, I deleted my answer.

[Caitlin Bemis]

Thank you, anyways!

Answer 1

[Bruno Luong]

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))

Number of intersections = 135

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

[Caitlin Bemis]

Thank you! This worked perfectly

[Caitlin Bemis]

Is there a way to see where the index is for each spot? That way I can segment the data into sections.

[Bruno Luong]

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]

Okay. I am not sure how to add things to the array. I attempted to add them to the seg2poly, but the response is either "too many input arguments" or "won't be used'. Unforunately, I am very new to Matlab trying something much more complex than my skills, obviously.

[Bruno Luong]

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