Phase shift of two sine curves

62 次查看（过去 30 天）
Bruce Rogers
Bruce Rogers 2021-6-30
评论： Mathieu NOE 2021-7-1
Hello everyone,
im trying to get the phase shift of two sine curves, one is the perfect sine wave (y=sine(x)) and the other one describes a robot movement. The robot movement is sometimes faster, but often slower than the theoretical sine wave. I'm trying to figure out how to get the phase shift, but I often get the error message "Index exceeds the number of array elements (10)" because the robot movement gives me sometimes more or sometimes less then 10 y-crossings. I'm comparing the nearest numbers to the y Crossing and I think the problem is, that Matlab thinks more numbers are the y Crossing although there only can be 10. Is there a function I can use?
I tried to split the data vector into 10 parts and search for the y-Crossing with mink funktion, but then also I get two points, that are really close to each other. Here a part of my code, the robot movement data is attached. It can be assumed that the first robot data value is the first y-Crossing, it doesn't matter which number the first data value has.
z_irl_temp = z_irl;
anzahl_z_irl = numel(z_irl_temp);
teilbar = round(anzahl_z_irl * 0.1); %Number of elements in vector * 0.1 so I get a Matrix with 10 columns
amnachsten = anzahl_z_irl + (teilbar - rem(anzahl_z_irl,teilbar)); %This is so I can get the number of cevtor rows that are divisible
for i = numel(z_irl):amnachsten
z_irl_temp(i) = 0; %filling the missung elements with zeros
end
z_irl_temp_numel = numel(z_irl_temp);
z_irl_temp_numel = z_irl_temp_numel/10;
neu = reshape(z_irl_temp,teilbar,[]); %reshape the new matrix to a matrix with 10 columns
const_irl = zeros(1,length(teilbar));
for i = 1:length(teilbar)
const_irl(i) = 389.387; %This is where the sine wave crosses the x-axis
end
const = const';
for i = 1:width(neu)
[~,idx_irl(:,i)] = mink(abs(neu(:,i)-const_irl),1);
end
for i = 1:width(neu)
if i > 1
idx_irl_temp(i) = idx_irl(i) + ((i-2) * teilbar);
time_yCrossings(i) = idx_irl_temp(i) * ti;
%idx_
else
time_yCrossings(i) = idx_irl(i) * ti;
end
end
Here's a picture of what exactly I'm trying to find. The green line describes the ideal sine curve, the blue one the real robot data.
I'm using Matlab R2020b for academic use. Thanks for your support, I'm thankful for your help!

请先登录，再回答此问题。

采纳的回答

Mathieu NOE
Mathieu NOE 2021-6-30
hello Bruce
the code you posted is lacking probably the section that load the data ; I could not really reproduce your plot.
so I prefered to show you my tools on dummy sinewaves - you can easily adapt it to your data.
this code generates two slighttly different sinewaves and compute the zero crossing time values for both positive and negative slope cases. Then it's fairly easy to compute the time delta between the two sets of points are shown at the end :
clc
clearvars
% dummy data
n=100;
x=(0:n-1)/n;
y1 = sin(30*x -0.5)+ 0.01*randn(1,n);
y2 = sin(31*x -0.7)+ 0.01*randn(1,n);
threshold = 0; % your value here
[t0_pos1,s0_pos1,t0_neg1,s0_neg1]= crossing_V7(y1,x,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points
[t0_pos2,s0_pos2,t0_neg2,s0_neg2]= crossing_V7(y2,x,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points
% ind => time index (samples)
% t0 => corresponding time (x) values
% s0 => corresponding function (y) values , obviously they must be equal to "threshold"
figure(1)
plot(x,y1,'b',t0_pos1,s0_pos1,'db',t0_neg1,s0_neg1,'+b','linewidth',2,'markersize',12);grid on
hold on
plot(x,y2,'g',t0_pos2,s0_pos2,'dg',t0_neg2,s0_neg2,'+g','linewidth',2,'markersize',12);grid on
hold off
legend('signal 1','signal 1 positive slope crossing points','signal 1 negative slope crossing points',...
'signal 2','signal 2 positive slope crossing points','signal 2 negative slope crossing points' );
% time difference for positive slope crossing points
dt_pos_slope = t0_pos1 - t0_pos2;
% time difference for negative slope crossing points
dt_neg_slope = t0_neg1 - t0_neg2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [t0_pos,s0_pos,t0_neg,s0_neg] = crossing_V7(S,t,level,imeth)
% [ind,t0,s0,t0close,s0close] = crossing_V6(S,t,level,imeth,slope_sign) % older format
% CROSSING find the crossings of a given level of a signal
% ind = CROSSING(S) returns an index vector ind, the signal
% S crosses zero at ind or at between ind and ind+1
% [ind,t0] = CROSSING(S,t) additionally returns a time
% vector t0 of the zero crossings of the signal S. The crossing
% times are linearly interpolated between the given times t
% [ind,t0] = CROSSING(S,t,level) returns the crossings of the
% given level instead of the zero crossings
% ind = CROSSING(S,[],level) as above but without time interpolation
% [ind,t0] = CROSSING(S,t,level,par) allows additional parameters
% par = {'none'|'linear'}.
% With interpolation turned off (par = 'none') this function always
% returns the value left of the zero (the data point thats nearest
% to the zero AND smaller than the zero crossing).
%
% [ind,t0,s0] = ... also returns the data vector corresponding to
% the t0 values.
%
% [ind,t0,s0,t0close,s0close] additionally returns the data points
% closest to a zero crossing in the arrays t0close and s0close.
%
% This version has been revised incorporating the good and valuable
% bugfixes given by users on Matlabcentral. Special thanks to
% Howard Fishman, Christian Rothleitner, Jonathan Kellogg, and
% Zach Lewis for their input.
% Steffen Brueckner, 2002-09-25
% Steffen Brueckner, 2007-08-27 revised version
% Copyright (c) Steffen Brueckner, 2002-2007
% brueckner@sbrs.net
% M Noe
% added positive or negative slope condition
% check the number of input arguments
error(nargchk(1,4,nargin));
% check the time vector input for consistency
if nargin < 2 | isempty(t)
% if no time vector is given, use the index vector as time
t = 1:length(S);
elseif length(t) ~= length(S)
% if S and t are not of the same length, throw an error
error('t and S must be of identical length!');
end
% check the level input
if nargin < 3
% set standard value 0, if level is not given
level = 0;
end
% check interpolation method input
if nargin < 4
imeth = 'linear';
end
% make row vectors
t = t(:)';
S = S(:)';
% always search for zeros. So if we want the crossing of
% any other threshold value "level", we subtract it from
% the values and search for zeros.
S = S - level;
% first look for exact zeros
ind0 = find( S == 0 );
% then look for zero crossings between data points
S1 = S(1:end-1) .* S(2:end);
ind1 = find( S1 < 0 );
% bring exact zeros and "in-between" zeros together
ind = sort([ind0 ind1]);
% and pick the associated time values
t0 = t(ind);
s0 = S(ind);
if strcmp(imeth,'linear')
% linear interpolation of crossing
for ii=1:length(t0)
%if abs(S(ind(ii))) > eps(S(ind(ii))) % MATLAB V7 et +
if abs(S(ind(ii))) > eps*abs(S(ind(ii))) % MATLAB V6 et - EPS * ABS(X)
% interpolate only when data point is not already zero
NUM = (t(ind(ii)+1) - t(ind(ii)));
DEN = (S(ind(ii)+1) - S(ind(ii)));
slope = NUM / DEN;
slope_sign(ii) = sign(slope);
t0(ii) = t0(ii) - S(ind(ii)) * slope;
s0(ii) = level;
end
end
end
% extract the positive slope crossing points
ind_pos = find(sign(slope_sign)>0);
t0_pos = t0(ind_pos);
s0_pos = s0(ind_pos);
% extract the negative slope crossing points
ind_neg = find(sign(slope_sign)<0);
t0_neg = t0(ind_neg);
s0_neg = s0(ind_neg);
end
  4 个评论
显示 2更早的评论
Bruce Rogers
Bruce Rogers 2021-7-1
Hello Mathieu,
I just used some easy if loop to check if the first crossing was detected by Matlab or not, so the whole code is running now \o/
Thanks a lot!

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Loops and Conditional Statements 的更多信息

标签

产品


版本

R2020b

Community Treasure Hunt

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

Start Hunting!

Translated by