How can I make the code make an action during part of the curve and make another action during certain part of the curve

1 view (last 30 days)
% if you take lock at the picture, I need to do an action durnig the main
% loop of the curve and need to do another action only during the minor
% loop 1,2 in the figure
load('m3rdr.mat') % this is the data of the picture
x=m3rdr(:,1);y=m3rdr(:,2);
plot(x,y)

Accepted Answer

Mathieu NOE
Mathieu NOE on 16 Mar 2022
hello Hassan
this is my suggestion - was not straightforward but hopefully it fills your need !
the minor loops data are xx and yy (cell arrays as they are not of same length)
clearvars
clc
% if you take lock at the picture, I need to do an action durnig the main
% loop of the curve and need to do another action only during the minor
% loop 1,2 in the figure
load('m3rdr.mat') % this is the data of the picture
x=m3rdr(:,1);y=m3rdr(:,2);
% computation of the curvature - need a bit of smoothing first
xs = smoothdata(x,'gaussian',15);
ys = smoothdata(y,'gaussian',15);
samples = 1:numel(x);
%% computation of the curvature
[curv] = compute_curvature(xs,ys);
[val,ind] = findpeaks(curv,'NPeaks',4,'SortStr','descend'); % find the 4 biggest peaks of curvature data
ind = sort(ind); % sort them
%% plots
figure(1)
plot(x,y);hold on
shrink_points = 2; % remove two leading and trailing points of the selected data
for ci =1:numel(ind)/2
ii = (1:2)+2*(ci-1)
ind2 = ind(ii);
dd = abs(ind2(2)-ind2(1)); % add the indexes of second half of data (based on ind2)
xx{ci} = x(ind2(1)+shrink_points:ind2(2)+dd-shrink_points);
yy{ci} = y(ind2(1)+shrink_points:ind2(2)+dd-shrink_points);
plot(xx{ci},yy{ci},'*')
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [curvature] = compute_curvature(x,y)
% run along the curve and find the radius of curvature at each location.
numberOfPoints = length(x);
curvature = zeros(1, numberOfPoints);
for t = 1 : numberOfPoints
if t == 1
index1 = numberOfPoints;
index2 = t;
index3 = t + 1;
elseif t >= numberOfPoints
index1 = t-1;
index2 = t;
index3 = 1;
else
index1 = t-1;
index2 = t;
index3 = t + 1;
end
% Get the 3 points.
x1 = x(index1);
y1 = y(index1);
x2 = x(index2);
y2 = y(index2);
x3 = x(index3);
y3 = y(index3);
a = sqrt((x1-x2)^2+(y1-y2)^2); % The three sides
b = sqrt((x2-x3)^2+(y2-y3)^2);
c = sqrt((x3-x1)^2+(y3-y1)^2);
A(t) = 1/2*abs((x1-x2)*(y3-y2)-(y1-y2)*(x3-x2)); % Area of triangle
curvature(t) = 4*A(t)/(a*b*c); % Curvature of circumscribing circle
end
end

More Answers (0)

Community Treasure Hunt

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

Start Hunting!