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

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

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

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,'.')

Figure contains an axes. The axes contains an object of type line.

Чтобы повторно отобразить неоднородно дискретизированный сигнал, можно вызвать 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])

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled.

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

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

Алгоритм преобразования в 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])

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled.

Отсутствующий сегмент соединяется линейной интерполяцией. Линейная интерполяция является методом по умолчанию, используемым 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])

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled using 'spline'.

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

По умолчанию 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,'.-')

Figure contains an axes. The axes contains an object of type line.

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

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

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

Figure contains an axes. The axes contains 2 objects of type line. These objects represent 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)')

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled (custom P and Q).

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

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

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

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

Figure contains an axes. The axes contains an object of type line.

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

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

Figure contains an axes. The axes contains an object of type line.

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

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])

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled.

Это работает достаточно хорошо. Однако повторно дискретизированный сигнал может быть сглажен дополнительно путем обеспечения 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])

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled (custom filter).

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

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

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

Figure contains an axes. The axes contains 2 objects of type line. These objects represent original, resampled (custom filter).

Эти программные продукты возникают из-за 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)

Figure contains an axes. The axes contains an object of type line.

Детрендированный сигнал теперь имеет обе свои конечные точки около нуля, что уменьшает введенные переходные процессы. Функции 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')

Figure contains an axes. The axes contains 2 objects of type line. These objects represent Original, Resampled (detrended, custom filter).

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

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

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

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

См. также