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)

Categories

Find more on Parallel Computing in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!