Повторная дискретизация неоднородно дискретизированных сигналов к желаемой скорости

The resample функция позволяет вам преобразовать неоднородно дискретизированный сигнал в новую равномерную скорость.

Создать синусоидальную частоту 500 Гц, отобранную нерегулярно примерно на 48 кГц. Мы симулируем неправильность, добавляя случайные значения к равномерному вектору.

rng default
nominalFs = 48000;
f = 500;
Tx = 0:1/nominalFs:0.01;
irregTx = sort(Tx + 1e-4*rand(size(Tx)));
x = sin(2*pi*f*irregTx);
plot(irregTx,x,'.')

Чтобы повторно отобразить неоднородно дискретизированный сигнал, можно вызвать resample с временным вектором входа.

Следующий пример преобразует наш исходный сигнал в равномерную частоту 44,1 кГц.

desiredFs = 44100;
[y, Ty] = resample(x,irregTx,desiredFs);
plot(irregTx,x,'.-',Ty,y,'o-')
legend('Original','Resampled')
ylim([-1.2 1.2])

Вы можете увидеть, что наш повторно дискретизированный сигнал имеет ту же форму и размер, что и исходный сигнал.

Выбор метода интерполяции

Алгоритм преобразования в resample лучше всего работает, когда входные выборки максимально близки к регулярным интервалам, поэтому полезно наблюдать, что может произойти, когда из выборочных данных отсутствует раздел входных выборок.

Следующий пример вырезает второй гребень входа синусоиды и применяет повторную дискретизацию.

irregTx(105:130) = [];
x = sin(2*pi*f*irregTx);
[y, Ty] = resample(x,irregTx,desiredFs);

plot(irregTx,x,'. ')
hold on
plot(Ty,y,'.-')
hold off
legend('Original','Resampled')
ylim([-1.2 1.2])

Отсутствующий сегмент соединяется линейной интерполяцией. Линейная интерполяция является методом по умолчанию, используемым resample функция для повторной выборки неоднородно выборочных данных.

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

Для малошумных сигналов с низкой пропускной способностью сплайны могут быть очень эффективными при использовании для восстановления исходного сигнала. Чтобы использовать кубический сплайн во время повторной дискретизации, укажите метод интерполяции 'spline':

[y, Ty] = resample(x,irregTx,desiredFs,'spline');

plot(irregTx,x,'. ')
hold on
plot(Ty,y,'.-')
hold off
legend('Original','Resampled using ''spline''')
ylim([-1.2 1.2])

Управление сеткой интерполяции

По умолчанию resample создает промежуточную сетку, которая является близким рациональным приближением отношения между желаемой частотой дискретизации и средней частотой дискретизации сигнала.

Если раздел ваших входных выборок содержит высокочастотные компоненты, можно управлять интервалом промежуточной сетки, выбирая целочисленные коэффициенты, p и q, чтобы выбрать это рациональное соотношение.

Исследуйте переходную характеристику заниженного фильтра второго порядка, который колеблется со скоростью около 3 Гц:

w = 2*pi*3;
d = .1002;
z = sin(d);
a = cos(d); 

t = [0:0.05:2 3:8];

x = 1 - exp(-z*w*t).*cos(w*a*t-d)/a;
plot(t,x,'.-')

Дискретизация переходной характеристики происходит с высокой скоростью, где она колеблется, и с низкой скоростью, где она отсутствует.

Теперь повторно отобразите сигнал на частоте 100 Гц как раз с помощью настроек по умолчанию:

Fs = 100;
[y, Ty] = resample(x,t,Fs);
plot(t,x,'. ')
hold on
plot(Ty,y)
hold off
legend('Original','Resampled (default settings)')

Огибающая колебания в начале формы волны ослабляется и колеблется медленнее, чем исходный сигнал.

resampleпо умолчанию интерполируется в сетку с регулярно разнесенными интервалами, которые соответствуют средней скорости дискретизации входного сигнала.

avgFs = (numel(t)-1) /(t(end)-t(1))

avgFs = 5.7500

Частота дискретизации сетки должна быть выше, чем в два раза больше самой большой частоты, которую вы хотите измерить. Частота дискретизации сетки, 5,75 выборки в секунду, ниже частоты выборки Найквиста, 6 Гц, частоты вызывного сигнала.

Чтобы сделать сетку более высокой частотой дискретизации, можно задать целочисленные параметры, p и q. resample настраивает частоту дискретизации сетки на Q * Fs/P, интерполирует сигнал, а затем применяет его внутренний преобразователь частоты дискретизации (увеличение дискретизации на P и понижение дискретизации на Q), чтобы восстановить желаемую частоту дискретизации Fs. Использование rat для выбора p и q.

Поскольку звон колебания составляет 3 Гц, задайте сетку со скоростью дискретизации 7 Гц, что немного выше скорости Найквиста. Заголовок 1 Гц учитывает дополнительное содержимое частоты из-за разрушающейся экспоненциальной огибающей.

Fgrid = 7;
[p,q] = rat(Fs/Fgrid)

p = 100

q = 7

[y, Ty] = resample(x,t,Fs,p,q);
plot(t,x,'.')
hold on
plot(Ty,y)
hold off
legend('Original','Resampled (custom P and Q)')

Определение фильтра сглаживания

В следующем примере можно просмотреть выход дигитайзера, который измеряет настройку дросселя на двигателе самолета. Установка дросселя неоднородно дискретизируется с номинальной скоростью 100 Гц. Мы постараемся повторно отобразить этот сигнал с равномерной частотой 10 Гц.

Вот выборки нашего исходного сигнала.

load engineRPM
plot(t,x,'.')
xlabel('Time (s)')
ylabel('RPM')

Наш сигнал квантован. Теперь увеличьте масштаб возрастающей области за временной интервал с 20 секунд до 23 секунд:

plot(t,x,'.')
xlim([20 23])

Сигнал медленно изменяется в этой области. Это позволяет вам удалить некоторый шум квантования с помощью фильтра сглаживания в ресамплере.

desiredFs = 10;
[y,ty] = resample(x,t,desiredFs);

plot(t,x,'.')
hold on
plot(ty,y,'.-')
hold off
xlabel('Time')
ylabel('RPM')
legend('Original','Resampled')
xlim([20 23])

Это работает достаточно хорошо. Однако повторно дискретизированный сигнал может быть сглажен дополнительно путем обеспечения resample фильтр с низкой частотой отключения.

Во-первых, установите интервал между сетками, чтобы быть примерно нашей номинальной частотой дискретизации 100 Гц.

nominalFs = 100;

Далее определите разумное p и q для получения желаемой скорости. Поскольку номинальная скорость 100 Гц и наша желаемая скорость 10 Гц, нужно децимировать на 10. Это эквивалентно настройке p значение 1 и настройка q по 10.

p = 1;
q = 10;

Можно поставить resample с вашим собственным фильтром. Чтобы иметь правильное временное выравнивание, фильтр должен быть нечетной длины. Длина фильтра должна быть в несколько раз больше p или q (в зависимости от того, какая величина больше). Установите частоту отключения равной 1/ q требуемая отсечка и затем умножение получившихся коэффициентов на p.

% ensure an odd length filter
n = 10*q+1;

% use .25 of Nyquist range of desired sample rate
cutoffRatio = .25;

% construct lowpass filter 
lpFilt = p * fir1(n, cutoffRatio * 1/q);

% resample and plot the response
[y,ty] = resample(x,t,desiredFs,p,q,lpFilt);

plot(t,x,'.')
hold on
plot(ty,y)
hold off
xlabel('time')
ylabel('RPM')
legend('Original','Resampled (custom filter)','Location','best')
xlim([20 23])

Удаление эффектов конечной точки

Теперь увеличьте изображение, чтобы просмотреть наш исходный сигнал. Обратите внимание, что в конечных точках имеется значительное смещение.

plot(t,x,'.',ty,y)
xlabel('time')
ylabel('RPM')
legend('original','resampled (custom filter)','Location','best')

Эти программные продукты возникают из-за resample принимает, что сигнал равен нулю вне границ сигнала. Чтобы уменьшить эффект этих разрывов, вычтите линию между конечными точками сигнала, выполните повторную дискретизацию, а затем добавьте назад линию к исходной функции. Вы можете сделать это, вычислив наклон и смещение линии между первой и последней выборкой, и используя polyval чтобы создать линию для вычитания.

% compute slope and offset (y = a1 x + a2)
a(1) = (x(end)-x(1)) / (t(end)-t(1));
a(2) = x(1);

% detrend the signal
xdetrend = x - polyval(a,t);
plot(t,xdetrend)

Детрендированный сигнал теперь имеет обе свои конечные точки около нуля, что уменьшает введенные переходные процессы. Функции resample а затем добавить назад тренд.

[ydetrend,ty] = resample(xdetrend,t,desiredFs,p,q,lpFilt);

y = ydetrend + polyval(a,ty);

plot(t,x,'.',ty,y)
xlabel('Time')
ylabel('RPM')
legend('Original','Resampled (detrended, custom filter)','Location','best')

Сводные данные

В этом примере показано, как использовать resample преобразование равномерно и неоднородно выбранные сигналы в фиксированную скорость.

Дальнейшее чтение

Для получения дополнительной информации о восстановлении неоднородно разнесённых выборок с помощью пользовательских сплайнов, можно ознакомиться с документацией Curve Fitting Toolbox™.