extracting one cycle from multiple sinusoidal waves

I am trying to extract one single cycle from each of the four responses attached into the excel sheet. Eventually I want to combine those extracted cycles into one array, for example A will have four columns:

column 1: one cycle of response 1

column 2: one cycle of response 2

column 3: one cycle of response 3

column 4: one cycle of response 4

Thanks in advance,

Answer 1

[Mathieu NOE]

hello

I used zero crossing interpolation method to select the last cycle of the data

the extracted data are interpolated because the start and stop time values are themselves interpolated data for maximum time precision and therefore may not belong to the actual data

clc
clearvars

T = readtable('responses.xlsx');

t = T.time_sec_;
R1 = T.resp1;
R2 = T.resp2;
R3 = T.resp3;
R4 = T.resp4;

threshold = 0; % your value here

[t0_pos1,s0_pos1,t0_neg1,s0_neg1]= crossing_V7(R1,t,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points 
figure(1),plot(t,R1,t0_pos1,s0_pos1,'db',t0_neg1,s0_neg1,'dg','linewidth',2,'markersize',12);grid on
legend('signal 1','signal 1 positive slope crossing points','signal 1 negative slope crossing points'  );
% NB  : t0_pos1 is interpolated so maybe best solution is to interpolate
% the data on new time axis; last cycle of data is used
new_t1 = linspace(t0_pos1(end-1),t0_pos1(end),100);
new_R1 = interp1(t,R1,new_t1);
figure(2),plot(new_t1,new_R1);

% same logic on R2
[t0_pos2,s0_pos2,t0_neg2,s0_neg2]= crossing_V7(R2,t,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points 
figure(3),plot(t,R2,t0_pos2,s0_pos2,'db',t0_neg2,s0_neg2,'dg','linewidth',2,'markersize',12);grid on
legend('signal 2','signal 2 positive slope crossing points','signal 2 negative slope crossing points'  );
new_t2 = linspace(t0_pos2(end-1),t0_pos2(end),100);
new_R2 = interp1(t,R2,new_t2);
figure(4),plot(new_t2,new_R2);

% same logic on R3
[t0_pos3,s0_pos3,t0_neg3,s0_neg3]= crossing_V7(R3,t,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points 
figure(5),plot(t,R3,t0_pos3,s0_pos3,'db',t0_neg3,s0_neg3,'dg','linewidth',2,'markersize',12);grid on
legend('signal 3','signal 3 positive slope crossing points','signal 3 negative slope crossing points'  );
new_t3 = linspace(t0_pos3(end-1),t0_pos3(end),100);
new_R3 = interp1(t,R3,new_t3);
figure(6),plot(new_t3,new_R3);

% same logic on R4
[t0_pos4,s0_pos4,t0_neg4,s0_neg4]= crossing_V7(R4,t,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points 
figure(7),plot(t,R4,t0_pos4,s0_pos4,'db',t0_neg4,s0_neg4,'dg','linewidth',2,'markersize',12);grid on
legend('signal 4','signal 4 positive slope crossing points','signal 4 negative slope crossing points'  );
new_t4 = linspace(t0_pos4(end-1),t0_pos4(end),100);
new_R4 = interp1(t,R4,new_t4);
figure(8),plot(new_t4,new_R4);

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


% 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 ~isempty(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);

else

    % empty output
    ind_pos = [];
    t0_pos = [];
    s0_pos = [];

    % extract the negative slope crossing points 
    ind_neg = [];
    t0_neg = [];
    s0_neg = [];
    
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % Addition:
% % Some people like to get the data points closest to the zero crossing,
% % so we return these as well
% [CC,II] = min(abs([S(ind-1) ; S(ind) ; S(ind+1)]),[],1); 
% ind2 = ind + (II-2); %update indices 
% 
% t0close = t(ind2);
% s0close = S(ind2);
end

[Abdulla Almesmari]

Thanks Mathieu This was perfect. In fact I have around 2000 responses and I believe it will be a time cosuming process to extract the single cycles for each one of them. Do you advise a better way?

[Mathieu NOE]

hello again

we can make a for loop code of it , depends now if your 2000 data comes from one file or 2000 files (format ? )

I can further help you in this task if you can share a few files for testing purpose

[Abdulla Almesmari]

In the workspace output, I have 1625 data obtained from simulink. I want to extract both columns from the data set (uy) and obtain only a single cycle for uy(:,1) and uy(:,2).

Is there a time-efficeint way to save those 1625 single-cycles in a separate vecotr?

I could not attach the file since it exceeds the max. size, please find the output file here:

https://we.tl/t-iDXwXgUgfW

[Mathieu NOE]

did you export the full simulation time in this file ?

maybe you can already reduce the file size by using the options of exporting only the last 1000 points or you shrink the file from matlab itself

[Mathieu NOE]

sorry, I meant even 1 second of data (the last second) would suffice

as your sampling frequency is 100 hz, this means the array length = 100 is enough to do the job

caan you see if you can save only the last second data ?

all the best

[Abdulla Almesmari]

sorry about that I only have the output data as I am attaching here:

https://we.tl/t-iDXwXgUgfW

The output from simulink were exported into a such format

[Mathieu NOE]

hello again

this is the code for your mat file

I used the first signal (square wave) to make the zero crossing detection action and select the last cycle

I exported the last cycle with time and both channels data. if you want only the data you can easily modify the code . rememnber that I interpolated the last cycle data on 100 samples , so all experiment data have same dimensions but may have different sampling times now

all experiments last cycles are concatenated in one array (data_out)

the other option would have been to select the actual data between the time index defined by the zero crossing x values , then we keep the original ( same ) sampling time but every single cycle data will have different length;

tell me what solution you prefer

clc
clearvars

load('matlab.mat');

[m,simulations] = size(out);

zcr_threshold = 0; % zero crossin detection threshold / your value here
data_out = [];

for ci = 1:simulations
    current_sim = out(1,ci);
    time = current_sim.tout;
    dt = mean(diff(time));
    Fs = 1/dt;
    data1 = current_sim.uy(:,1); %  I want to extract both columns from the data set (uy)  
    data2 = current_sim.uy(:,2); %  and obtain only a single cycle for uy(:,1) and uy(:,2). 
    
 
    % do zero crossing detection on first signal (square wave) and extract last data cycle  
    [t0_pos1,s0_pos1,t0_neg1,s0_neg1]= crossing_V7(data1,time,zcr_threshold,'linear'); % positive (pos) and negative (neg) slope crossing points 
    if ci<11 % plot only the 10 first experiments data
        figure(ci),subplot(211),plot(time,data1,time,data2,t0_pos1,s0_pos1,'db',t0_neg1,s0_neg1,'dg','linewidth',2,'markersize',12);grid on
        legend('signal 1','signal 2','signal 1 positive slope crossing points','signal 1 negative slope crossing points'  );
    end
    % NB  : t0_pos1 is interpolated so maybe best solution is to interpolate
    % the data on new time axis; last cycle of data is stored in array
    % "data_out"
    new_t1 = linspace(t0_pos1(end-1),t0_pos1(end),100);
    new_data1 = interp1(time,data1,new_t1);
    new_data2 = interp1(time,data2,new_t1);
    subplot(212),plot(new_t1,new_data1,new_t1,new_data2);
    
    tmp = [new_t1(:) new_data1(:) new_data2(:)]; % one experiment data format : time / data 1 / data 2
    data_out = [data_out tmp];      % vertical concatenation of all experiments data 
    
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).
%


% 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 ~isempty(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);

else

    % empty output
    ind_pos = [];
    t0_pos = [];
    s0_pos = [];

    % extract the negative slope crossing points 
    ind_neg = [];
    t0_neg = [];
    s0_neg = [];
    
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % Addition:
% % Some people like to get the data points closest to the zero crossing,
% % so we return these as well
% [CC,II] = min(abs([S(ind-1) ; S(ind) ; S(ind+1)]),[],1); 
% ind2 = ind + (II-2); %update indices 
% 
% t0close = t(ind2);
% s0close = S(ind2);
end

[Abdulla Almesmari]

Thanks Mathieu your solution works perfectly

[Mathieu NOE]

My pleasure !