Подсчет сегментов имульсов в сигнале MATLAB

Уважаемые коллеги, добрый вечер! В общем, возникла проблема следующего характера.

Имеется сигнал, достаточно большой объем точек, длительность порядка 35-40 секунд. Он представлят собой последовательность импульсов тока, длительностью примерно 5 мкс. Есть, например, такой сигнал, где последовательности сигналов (импульсов) составляют такие участки (sig1.mat):

         

Суть задачи такова, что необходимо подсчитывать количество таких сегментов (выделены рамкой), и количесто импульсов в кажом таком сегменте, собственно говоря, просуммировав весь массив получим общее число рабочих импульсов. Для начала, я хотел использовать конструкцию счетчика путем применения логической индекасации, полагая, что однозначно, импульс по амплитуде имеет значение больше 200. Тогда можно для каждого импульса брать, определять начальный пик имульса, затем находить разницу до следующего пика, который больше порогового значения (200) и как только дальше идет промежуток, когда значение около нуля, т.е интервал следования между импульсами (больше нуля) использовать nnz. Т.е что-то примерно такое:

numberOfPulses = nnz(diff(outSig_1_curr > 1) > 0);

Однако, в сигнале могут присутствовать и такие фрагменты файл sig16.mat:

Тогда попрбовал напрямую использовать индексацию, и вычислять интервалы вручную (к сожалению в этом случае минимальный размер окна, т.е. длительность такого фрагмента 0.1. мс и расстояние между импульсами 0.1 мс):

clc
clearvars
load('sig1.mat') % outSig_1_curr2
samples = numel(outSig_1_curr2);
Fs = 1000000; %  Hz
dt = 1/Fs;
t = (0:samples-1)*dt;
% % extract signal between t = 18 and t = 20 s
% ind1 = round(0.001*Fs);
% ind2 = round(0.1*Fs);
% t = t(ind1:ind2);
% outSig_1_curr2 = outSig_1_curr2(ind1:ind2);
% use tmp instead of outSig_1_curr2 data for islocalmax 
% tmp is the same as outSig_1_curr2 , except that negative values are
% replaced by zeroes
tmp = outSig_1_curr2;
tmp(tmp<0) = 0;
%% find peaks
% find positive local maxima
[tf, P] = islocalmax(tmp,'MinProminence',200); % see remark above
t_peak = t(tf);
outSig_1_peak = outSig_1_curr2(tf);
%% localize / extract valid data blocks 
blocks_distance_samples = 100; % select data blocks  that are distant above value = blocks_distance_samples (in samples)
min_contiguous_samples = 100; % select data blocks  that are contiguous with min amount of samples = min_contiguous_samples
detect_signal = zeros(size(t))+0.5; % initialisation
    % now define start en end point of  segments above threshold
        %%%%%%%%%%
        % This locates the beginning /ending points of data groups
        ind = find(tf); % logical to index conversion
        D = diff([0;ind;0]);
        ind_ends = find(D > blocks_distance_samples) - 1;
        ind_ends = ind_ends(ind_ends>0); % remove potential zero value
        ends = ind(ind_ends);
        begin = ind(ind_ends+1);
        % remove incorrect begin / ends points (like end point before begin  point and vice versa)
        ends = ends(ends>begin(1));
        begin = begin(begin<ends(end));
        m = min(numel(begin),numel(ends));
        ends = ends(1:m);
        begin = begin(1:m);
        
        
        %%%%%%%%%%
        
    length_ind = ends - begin; % compute block length
    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 
    
    t_peak2 = t(ind);
    outSig_1_peak2 = outSig_1_curr2(ind);
    % define the begin / ending x, y values of raw data 
    t_peak2_begin = t(begin);
    outSig_1_curr2_begin = outSig_1_curr2(begin);
    t_peak2_ends = t(ends);
    outSig_1_curr2_ends = outSig_1_curr2(ends);
    
for ci = 1:length(begin)
    ind = (t>=t_peak2_begin(ci) & t<=t_peak2_ends(ci));
    detect_signal(ind) = 200;
    % include here your code if you want to extract the data blocks
    % individually
    ind4 = (t_peak>=t_peak2_begin(ci) & t_peak<=t_peak2_ends(ci));
    number_of_peaks(ci) = numel(find(ind4));
    
end
number_of_peaks
% plot
figure(1)
plot(t,outSig_1_curr2,t_peak,outSig_1_peak,'dr',t(begin),outSig_1_curr2(begin),'dg',t(ends),outSig_1_curr2(ends),'dk',t,detect_signal,'k');
legend('signal',' peaks ','begin points','end points','selected signal')

Однако, к сожалению, в этом случае очень тяжело распознать фргаменты со второго рисунка и подсчитывать такие фрагменты. Между импульсами в таких фрагментах длительность почти нулевая, и задать размер между импульсами в 1-2 отсчета не получается, происходит ошибка индексации.

Каким образом, можно задать расстояние в 2-3 сэмпла и при этом различать такие фрагменты? Или можно же использовать конструкцию с nnz, и в полученном массиве логических единиц определять общее число импульсов и имульсов фрагмента?

К примеру:

idx1 = [0 0 0 1 1 1 0 0 1 1 1 0 0 0 0 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 0 0];
diffs = diff([0 idx1 0]);

% первый индекс
loc = find(diffs == 1);

Ns = find(diffs == -1) - loc;

C = mat2cell([loc' Ns'],ones(size(loc)))

Может быть есть готовые решения или более корректное использование функций islocalmax, islocalmin для подсчета ипульсов? Возможно ли задать параметр 'MinProminence как дииапазон (например, часть фргаментов, которые на втором рисунке лежит примерно выше 1500 по амлитуде, а каждый пик импульса такого фргамента приблизительно выше порогового значения max(signal)/1.03321). На данный момент я не могу предположить каким образом, можно подсчитать количесвто таких импульсов в каждом сегменте и количество таких сегментов.

Сами сигналы доступны по ссылке:https://drive.google.com/drive/folders/1NB3QQXVV023fa4fdK6MzHqLHI0WH3aLC?usp=sharing

Отмечу, что выделены черным цветом - это одна группа импульсов, рабочие, а те, что красным - это другая гурппа сегментнов (короткое замыкание)

                   

Однако, я получаю все равно ошибку индексации:

Array indices must be positive integers or logical values.

Error in movavg (line 43)
        ends = ind(ind_ends);

freq = 1e6;
step = 1/freq;
time_sig = sig16(:,1);
data = sig16(:,2);
 
fl1 = 120;
[up1,lo1] = envelope(data,1000,'analytic');
upp1 = filter(ones(1,fl1)/fl1,1,up1);
upp1_norm = circshift(upp1,-95);
 
 
% upp1_norm - огибающая
[tf, P] = islocalmin(upp1_norm,'MinProminence',min(upp1_norm));
t_peak = time_sig(tf);
outSig_1_peak = upp1_norm(tf);
 
 
%% localize / extract valid data blocks
blocks_distance_samples = 100; % длина блока, состоящие из импульсов, где огибающая ниже средней линии
min_contiguous_samples = 100; % рамзер окна
detect_signal = zeros(size(time_sig))+0.5;
 
    % старт и конец блока с определенным пороговым значением огибающей min(upp1_norm)
   
        %%%%%%%%%%
        % This locates the beginning /ending points of data groups
        ind = find(tf);
        D = diff([0;ind;0]);
        ind_ends = find(D > blocks_distance_samples) - 1;
        ends = ind(ind_ends);
        begin = ind(ind_ends+1);
        % удалить начальные и конечне некорректные точки сегмента
        ends = ends(ends>begin(1));
        begin = begin(begin<ends(end));
        m = min(numel(begin),numel(ends));
        ends = ends(1:m);
        begin = begin(1:m);
       
       
        %%%%%%%%%%
       
    length_ind = ends - begin; % вычисляем длину блока
    ind2= length_ind>min_contiguous_samples; % проверка корректности длины
    begin = begin(ind2); % выбранные точки
    ends = ends(ind2); % выбранные точки
   
    t_peak2 = t(ind);
    outSig_1_peak2 = upp1_norm(ind);
 
    t_peak2_begin = t(begin);
    outSig_1_curr2_begin = upp1_norm(begin);
    t_peak2_ends = t(ends);
    outSig_1_curr2_ends = upp1_norm(ends);
   
for ci = 1:length(begin)
    ind = (t>=t_peak2_begin(ci) & t<=t_peak2_ends(ci));
    detect_signal(ind) = min(upp1_norm);
 
    ind4 = (t_peak>=t_peak2_begin(ci) & t_peak<=t_peak2_ends(ci));
    number_of_peaks(ci) = numel(find(ind4));
   
end
 
number_of_peaks;

Но пока, этот вариант работает не совсем хорошо...

Так что я буду благодарен за помощь! Огромное спасибо!