Заключите сигнал, получив верхнюю и нижнюю огибающую в MATLAB

У меня есть сигнал, который колеблется вокруг некоторых значений. Я хочу заключить сигнал сверху и снизу, как показано на следующем рисунке: example

Если сигнал только растет или уменьшается, отлично работает работающий скрипт maxima:

function [out] = rMax(X)
    Y=X;
    maximum=X(1);
    for k=1:length(X)
        if X(k)<=maximum
        Y(k)=maximum;
        else
            maximum=X(k);
        end
    end
    out=Y;
end

Однако, когда сигнал меняется, я больше не могу использовать этот метод. Есть ли способ добиться этого в MATLAB или Mathematica?


person Jared Lo    schedule 11.07.2016    source источник
comment
Для этой или любой другой функции я настоятельно рекомендую вам НЕ использовать переменную с именем max, так как она затеняет встроенную функцию Matlab и может вызвать всевозможные проблемы с отладкой, особенно когда это делается по незнанию.   -  person Hoki    schedule 12.07.2016
comment
да я знаю, спасибо   -  person Jared Lo    schedule 13.07.2016


Ответы (2)


Вы можете использовать movmax и movmin в сочетании с соответствующими размерами окна n1 > n2 в зависимости от ваших потребностей. Эти две функции были представлены в R2016a. Посмотрите в нижней части моего ответа на самостоятельно написанную замену для них, которая работала до R2016a.

Чтобы получить верхние границы, вы можете использовать:

xMax = movmax(x,[n1,n1]);
xMax = movmin(xMax,[n2,n2]);

Для нижних границ просто поменяйте местами movmin и movmax:

xMin = movmin(x,[n1,n1]);
xMin = movmax(xMin,[n2,n2]);

Результат может выглядеть так:

результат

Если вы выберете n2=n1, границы будут очень узкими, когда в данных есть пик x. Если вы выберете большую разницу, просто сделав n2 меньше, чем n1, вы получите более длинную прямую линию на пике. Побочным эффектом является то, что границы начинают отдаляться, когда сигнал переходит от 1 к -1 или наоборот. Если n1 выбрано слишком большим, небольшая часть сигнала около 100 не будет извлечена.


Вот полный код для создания приведенного выше рисунка с некоторыми примерами данных, соответствующими графику в вашем вопросе. Просто поэкспериментируйте со значениями n1 и n2, чтобы увидеть их эффект.

n1 = 20;        % for first window
n2 = 18;        % for second window

% generate sample data
t = 20:0.1:220;
x = -ones(size(t));
x(t>60&t<100) = 1;
x(t>105&t<135) = 1;
x(t>145&t<155) = 1;
x = x + 0.4*randn(size(x));

% get upper bounds
xMax = movmax(x,[n1,n1]);
xMax = movmin(xMax,[n2,n2]);

% get lower bounds
xMin = movmin(x,[n1,n1]);
xMin = movmax(xMin,[n2,n2]);

% draw figure for illustration
figure; hold on;
plot(t,x);
xlim([20,220]);
ylim([-3,3]);
plot(t,xMax,'r','LineWidth',1.1);
plot(t,xMin,'Color',[0,0.5,0],'LineWidth',1.1);

До R2016a

Чтобы иметь базовую функциональность movmin и movmax в версиях MATLAB до R2016a, мы можем реализовать наши собственные функции. Следовательно, нам нужно применить скользящий минимум (или скользящий максимум соответственно), который может быть достигнут довольно легко. Чтобы сохранить совместимость с функциями в R2016a, мы собираемся реализовать случай, когда второй аргумент является скаляром, а вектор — вектором из двух элементов. Это относится к следующему синтаксису, аналогичному показанному здесь, с тем ограничением, что x должен быть вектор.

y = movingmax(x,k)
y = movingmax(x,[kb kf])

Вот код для movingmax, который заменяет movmax:

function y = movingmax(x,n)
if length(n) > 1
    a=n(1); b=n(2);
else
    b=floor((n-1)/2); a=b+mod(n-1,2);
end
s = size(x);
xp = [ones(a,1)*x(1);x(:);ones(b,1)*x(end)];
y = zeros(size(x));
for k = 1:length(x)
    y(k) = max(xp(k:k+a+b));
end
y = reshape(y,s);

Вот код для movingmin, который заменяет movmin:

function y = movingmin(x,n)
if length(n) > 1
    a=n(1); b=n(2);
else
    b=floor((n-1)/2); a=b+mod(n-1,2);
end
s = size(x);
xp = [ones(a,1)*x(1);x(:);ones(b,1)*x(end)];
y = zeros(size(x));
for k = 1:length(x)
    y(k) = min(xp(k:k+a+b));
end
y = reshape(y,s);
person Matt    schedule 12.07.2016
comment
Ваш ответ, похоже, делает именно то, что я искал. Плохо только то, что функции movmax и movmin работают для matlab2016a, у меня сейчас доступ только к 2015 версии. - person Jared Lo; 13.07.2016
comment
@JaredLo Я только что обновил свой ответ, включив в него решение, которое работает для версий до R2016a. См. часть внизу, которая предлагает замену для двух функций. Результаты точно такие же. - person Matt; 13.07.2016

Я добавляю второе условие к моим работающим функциям максимума/минимума.

function [out] =rMin(X)

Y=X;

minimum=X(1);

for k=2:length(X)

    if X(k)>=minimum
        Y(k)=minimum;
    else
        minimum=X(k);
    end

    if X(k)>=(Y(k-1)+3)
        minimum=X(k);
        Y(k)=minimum; 
    end
end 

out=Y;

end

и я понимаю, что это все еще не работает в некоторых частях сигнала введите здесь описание изображения

person Jared Lo    schedule 12.07.2016