How to count active faults in a specific time interval

hello ,i want to detect faults. if the data ( residual) goes beyond a certain threshold and comes back , it is not a fault so it shows 1 but if the residual goes the threshold for consecutive number of times with a specific time-interval then it shows 1 .When there is no fault , it should show 1 . The fault will be plotted against time . How do i implement this?

[Jon]

As stated you "show 1" in every case fault or no fault. Do you mean when there is no fault it should "show" 0?

[Kwaku Junior]

my mistake , i meant show 0 there is no fault

Answer 1

[Mathieu NOE]

hello

try this code ; it select segments of data above a given threshold and longer than min_contiguous_samples samples;

% dummy data
samples=10000;
time= 3*(0:samples-1)/samples;
dt = mean(diff(time));
data = max(0,0.03*time+ sin(6*time.^2 -0.5));
fault_signal = zeros(size(time));

%% parameters
min_contiguous_samples = 400; % select segments only if they are at least this length => faulty signal = 1
threshold = 0.5; % 1 = max (100%)  of data value

%% main loop %%%%

    ind = (data>threshold*max(data)); % find data above threshold
    % now define start en end point of  segments above threshold
    [begin,ends] = find_start_end_group(ind);
    length_ind = ends - begin;
    ind2= length_ind>min_contiguous_samples; % check if their length is valid (above min_contiguous_samples value)
    begin = begin(ind2); % selected  points 
    ends = ends(ind2); % selected  points 
    
    time2 = time(ind);
    data2 = data(ind);
    % define the begin / ending x, y values of raw data 
    time2_begin = time(begin);
    data_begin = interp1(time,data,time2_begin);
    time2_ends = time(ends);
    data_ends = interp1(time,data,time2_ends);
    
for ci = 1:length(begin)
    ind = (time>=time2_begin(ci) & time<=time2_ends(ci));
    fault_signal(ind) = 1;
end

figure(1),
subplot(211),plot(time,data,'k',time2,data2,'.r',time2_begin,data_begin,'*c',time2_ends,data_ends,'*m','MarkerSize',12);
xlabel('time');
ylabel('amplitude');
legend('signal',['signal above ' num2str(threshold*100) ' % threshold'] ,'begin points','end points');
subplot(212),plot(time,fault_signal,'k');
xlabel('time');
ylabel('amplitude');
legend('fault signal');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [begin,ends] = find_start_end_group(ind)
% This locates the beginning /ending points of data groups
    D = diff([0,ind,0]);
    begin = find(D == 1);
    ends = find(D == -1) - 1;
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function x_rms = my_rms(x)
    x_rms = sqrt(mean(x.^2));
    
end

[Kwaku Junior]

thank you very much , i got this error when i run it ,"Function definitions are not permitted in this

context."

[Mathieu NOE]

maybe you should try save each function in a separate m file and remove from main code

also you don't need this one (a relique of another code) :

function x_rms = my_rms(x)

x_rms = sqrt(mean(x.^2));

end

[Mathieu NOE]

even simpler - I put back the code from function into the main section

works fine on my side

% dummy data
samples=10000;
time= 3*(0:samples-1)/samples;
dt = mean(diff(time));
data = max(0,0.03*time+ sin(6*time.^2 -0.5));
fault_signal = zeros(size(time));

%% parameters
min_contiguous_samples = 400; % select segments only if they are at least this length => faulty signal = 1
threshold = 0.5; % 1 = max (100%)  of data value

%% main loop %%%%

    ind = (data>threshold*max(data)); % find data above threshold
    % now define start en end point of  segments above threshold
        %%%%%%%%%%
        % This locates the beginning /ending points of data groups
        D = diff([0,ind,0]);
        begin = find(D == 1);
        ends = find(D == -1) - 1;
        %%%%%%%%%%
    length_ind = ends - begin;
    ind2= length_ind>min_contiguous_samples; % check if their length is valid (above min_contiguous_samples value)
    begin = begin(ind2); % selected  points 
    ends = ends(ind2); % selected  points 
    
    time2 = time(ind);
    data2 = data(ind);
    % define the begin / ending x, y values of raw data 
    time2_begin = time(begin);
    data_begin = interp1(time,data,time2_begin);
    time2_ends = time(ends);
    data_ends = interp1(time,data,time2_ends);
    
for ci = 1:length(begin)
    ind = (time>=time2_begin(ci) & time<=time2_ends(ci));
    fault_signal(ind) = 1;
end

figure(1),
subplot(211),plot(time,data,'k',time2,data2,'.r',time2_begin,data_begin,'*c',time2_ends,data_ends,'*m','MarkerSize',12);
xlabel('time');
ylabel('amplitude');
legend('signal',['signal above ' num2str(threshold*100) ' % threshold'] ,'begin points','end points');
subplot(212),plot(time,fault_signal,'k');
xlabel('time');
ylabel('amplitude');
legend('fault signal');

[Kwaku Junior]

it is works well but when i use my data it doesnt work this is the error ( Error using horzcat Dimensions of matrices being concatenated are not consistent.) , i dont know if its because my data is array of 1s and 0s . also it is a 50000x1 double

[Mathieu NOE]

can you show where the error occurs in your code ?

[Kwaku Junior]

-Error using horzcat Dimensions of matrices being concatenated are not consistent.

D = diff([0,ind,0]);

this is where the error occurs

also the samples being used , how do i use my array ( .mat file), the problems seem to come from there.

[Mathieu NOE]

ok found the issue , as your data is column oriented and my demo was using row oriented data

simply change this

D = diff([0,ind,0]);

into this

D = diff([0;ind;0]);

full code :

% dummy data
samples=10000;
time= 3*(0:samples-1)'/samples;
dt = mean(diff(time));
data = round(max(0,0.03*time+ sin(6*time.^2 -0.5)));
fault_signal = zeros(size(time));

%% parameters
min_contiguous_samples = 400; % select segments only if they are at least this length => faulty signal = 1
threshold = 0.5; % 1 = max (100%)  of data value

%% main loop %%%%

    ind = (data>threshold*max(data)); % find data above threshold
    % now define start en end point of  segments above threshold
        %%%%%%%%%%
        % This locates the beginning /ending points of data groups
        D = diff([0;ind;0]);
        begin = find(D == 1);
        ends = find(D == -1) - 1;
        %%%%%%%%%%
    length_ind = ends - begin;
    ind2= length_ind>min_contiguous_samples; % check if their length is valid (above min_contiguous_samples value)
    begin = begin(ind2); % selected  points 
    ends = ends(ind2); % selected  points 
    
    time2 = time(ind);
    data2 = data(ind);
    % define the begin / ending x, y values of raw data 
    time2_begin = time(begin);
    data_begin = interp1(time,data,time2_begin);
    time2_ends = time(ends);
    data_ends = interp1(time,data,time2_ends);
    
for ci = 1:length(begin)
    ind = (time>=time2_begin(ci) & time<=time2_ends(ci));
    fault_signal(ind) = 1;
end

figure(1),
subplot(211),plot(time,data,'k',time2,data2,'.r',time2_begin,data_begin,'*c',time2_ends,data_ends,'*m','MarkerSize',12);
xlabel('time');
ylabel('amplitude');
legend('signal',['signal above ' num2str(threshold*100) ' % threshold'] ,'begin points','end points');
subplot(212),plot(time,fault_signal,'k');
xlabel('time');
ylabel('amplitude');
legend('fault signal');

[Kwaku Junior]

it works now , thank you very much

[Mathieu NOE]

My pleasure !

Answer 2

[Jon]

Here is another approach (uses repelem trick from @Jan https://www.mathworks.com/matlabcentral/answers/382011-how-to-count-the-number-of-consecutive-identical-elements-in-both-the-directions-in-a-binary-vecto?s_tid=srchtitle)

% set up a small example
thresh = 0.5 % threshold for fault
minCount = 2; % min number of repeated values to be considered a fault
r = rand(1,10) % example vector of residuals

% make logical array with values of residual greater than threshold set
% true
f = r > thresh

% mark beginning and ends of consecutive runs of zeros and ones
d = [true,diff(f)~=0,true]

% find run lengths
n = diff(find(d))

% make vector with run length corresponding to each element
runLength = repelem(n,n);

% just keep run lengths where it is over threshold
runLength = runLength.*f

% assign logical indices for points that are considered faults
isFault = runLength >= minCount

% illustrative plots
idx = 1:numel(r);
plot(idx,r,idx(~isFault),r(~isFault),'og',idx(isFault),r(isFault),'*r')

[Kwaku Junior]

Thank you

'Undefined function 'repelem' for input arguments of type 'double'."

i got this error from your code

[Jon]

Strange, runs without error on my computer, also you can see here that repelem is a standard MATLAB function https://www.mathworks.com/help/matlab/ref/repelem.html

Did the code above throw the error or was it your adaptation of the code. If it was your modified code perhaps you mispelled the function name. You can also try typing

which repelem

on the command line. On my computer I get

>> which repelem

C:\Program Files\MATLAB\R2022a\toolbox\matlab\elmat\repelem.m

Which shows that it is part of the standard MATLAB install

[Kwaku Junior]

repelem isnt working because my matlab is 2013

[Jon]

Ahh and it wasn't introduced until 2015. The approach above can be done with a little modification just using diffs and cumsums. I would have to think about that a little more though. Is it still of interest?

Otherwise, is there some reason why you can't update your MATLAB, 2013 is almost 10 years old, a lot of advances since then.

[Mathieu NOE]

hello again

there are "home made" alternatives to repelem if you don't have it :

repelem substitute for older matlab version - (mathworks.com)

see in the code below (expanded)

% set up a small example
thresh = 0.5; % threshold for fault
minCount = 5; % min number of repeated values to be considered a fault
r = rand(1,100) % example vector of residuals
% make logical array with values of residual greater than threshold set
% true
f = r > thresh;
% mark beginning and ends of consecutive runs of zeros and ones
d = [true,diff(f)~=0,true];
% find run lengths
n = diff(find(d));
% make vector with run length corresponding to each element
% runLength = repelem(n,n);
runLength = my_repelem1(n,n); % alternative #1
% runLength = my_repelem2(n,n); % alternative #2

% just keep run lengths where it is over threshold
runLength = runLength.*f;
% assign logical indices for points that are considered faults
isFault = runLength >= minCount;
% illustrative plots
idx = 1:numel(r);
plot(idx,r,idx(~isFault),r(~isFault),'og',idx(isFault),r(isFault),'*r')

%%%%%%%%%%%%%%%%%%%%%% alternatives to repelem %%%%%%%%%%%%%%%%%%%%%
function out = my_repelem1(A,R)
    %   B = my_repelem1(A, R), returns an array with each element of A
    %   repeated according to R. 
    out = cell2mat(arrayfun(@(a,r)repmat(a,1,r),A,R,'uni',0));

end

%%%%%%%%%%%%%%%%%%%%%%
function out = my_repelem2(A,R)
    %   B = my_repelem1(A, R), returns an array with each element of A
    %   repeated according to R. 
    out = cell(1,length(R)) ;
    for i = 1:length(R)
        out{i} = ones(1,R(i))*A(i) ;
    end
    out = [out{:}];
end