Problems in for loop with cell array

Hello everyone,

I'm trying to get the phase shift of a normal sine wave with different sine waves that were produced by a robot. Every run of the experiment has a different period, so its y=sine(factor * x). I used Crossing_V7 function that was given to me a few days ago. When I used it with only one robot sine wave it worked fine, but now I want to use it for all five attempts with different parameters and I get the error message, that the length is different. Does anyone of you see the mistake? I'm pretty new to cell array calculations.

clear;

factor = [0.25, 0.5, 1.0, 2.0, 4.0];

z_irl =  cell(1,5);

irl = readtable("C:\Users\Bruce Rogers\Documents\RoboDK\halfpi\Roboterfahrten\RSI_Auswertung_halfpi\min_vel12_halfpi.txt");

z_irl_temp = irl{:,4};

z_irl{:,1} = z_irl_temp;

irl = readtable("C:\Users\Bruce Rogers\Documents\RoboDK\1pi\Roboterfahrten\RSI_Auswertung_1pi\min_vel12_1pi.txt");

z_irl_temp = irl{:,4};

z_irl{:,2} = z_irl_temp;

irl = readtable("C:\Users\Bruce Rogers\Documents\RoboDK\2pi\Roboterfahrten\RSI_Auswertung_2pi\min_vel12_2pi.txt");

z_irl_temp = irl{:,4};

z_irl{:,3} = z_irl_temp;

irl = readtable("C:\Users\Bruce Rogers\Documents\RoboDK\4pi\Roboterfahrten\RSI_Auswertung_4pi\min_vel12_4pi.txt");

z_irl_temp = irl{:,4};

z_irl{:,4} = z_irl_temp;

irl = readtable("C:\Users\Bruce Rogers\Documents\RoboDK\8pi\Roboterfahrten\RSI_Auswertung_8pi\min_vel12_8pi.txt");

z_irl_temp = irl{:,4};

z_irl{:,5} = z_irl_temp;

%%

for j = 1:length(factor)
    
z_irl = z_irl{1,j};
    
num_irl = numel(z_irl);

t = (1:num_irl)*0.004;

A = 50;

b = 10; %Durchgänge

ti = 0.004;

%x = zeros(length(z_irl));

x(:,j) = 0:ti:(1/factor(j))*10*pi;

y(:,j) = -A*sin(factor(j)*x) + 389.387; 

%y(j) = y(j)';

num_th = numel(x);

T = (1:length(x))*ti;

threshold = 389.387; % your value here
[t0_pos1(:,j),s0_pos1(:,j),t0_neg1(:,j),s0_neg1(:,j)]= crossing_V7(y(:,j),T(:,j),threshold,'linear'); % positive (pos) and negative (neg) slope crossing points %%Here I get the error message
[t0_pos2(:,j),s0_pos2(:,j),t0_neg2(:,j),s0_neg2(:,j)]= crossing_V7(z_irl(:,j),t(:,j),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(7)
% plot(T,y,'g',t0_pos1,s0_pos1,'+g',t0_neg1,s0_neg1,'+g','linewidth',2,'markersize',12);grid on
% hold on
% plot(t,z_irl,'b',t0_pos2,s0_pos2,'+b',t0_neg2,s0_neg2,'+b','linewidth',2,'markersize',12);grid on
% hold off
% title('Nullstellen beider Wellen, Grün: Theorie, Blau: Realität')
% xlabel('time');
% ylabel('distance');



% time difference for positive slope crossing points
dt_pos_slope(j) = t0_pos1(j) - t0_pos2(j);
% time difference for negative slope crossing points

if numel(t0_neg1(j)) ~= 5

    t0_neg1_temp(j) = 0.004;

    t0_neg1(j) = [t0_neg1_temp(j),t0_neg1(j)];

end

if numel(t0_neg2(j)) ~= 5

    t0_neg2_temp(j) = 0.004;

    t0_neg2(j) = [t0_neg2_temp(j),t0_neg2(j)];

end

dt_neg_slope(j) = t0_neg1(j) - t0_neg2(j);

end

%%

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
% [e-mail]
% 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!');                      %This is the error message I get
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

Answer 1

[David Hill]

Why not just use fit() function

x=0:.01:2*pi;
y=(5+.05*randn(size(x))).*sin((3+.05*randn(size(x))).*x);
f=fit(x',y','sin1');
period=f.b1;

[Bruce Rogers]

Hi David,

i think you misunderstood my problem. I don't want to fit two sine waves together, I want the phase shift. This works if I load every experiment with one parameter individually, but later I want to generate a contour plot with the different periods, different parameters and the different phase shifts. In the image below you can see what I mean.

[David Hill]

To obtain the amplitude, period, and shift of each sine wave, just run the above. Not sure how your are storing your data for each sine wave (x,y).

for k=1:numOfsineWaves
  f=fit(x{k}',y{k}','sin1');%assumed the sine waves are stored in cell arrays since the size varies
  m(k,:)=[f.a1,f.b1,f.c1];%stores amplitude,period,shift in matrix for all sine waves
end

[Bruce Rogers]

Now I get what you mean. Thanks David, it helped, but I used another method, because otherwise I would had to change too much on the program