불균일하게 샘플링된 신호 리샘플링하기

불균일하게 샘플링된 신호를 원하는 레이트로 리샘플링하기

resample 함수를 사용하면 불균일하게 샘플링된 신호를 균일한 새 레이트로 변환할 수 있습니다.

약 48kHz로 불규칙적으로 샘플링된 500Hz 정현파를 생성합니다. 난수 값을 등간격 벡터에 추가하여 불규칙성을 시뮬레이션합니다.

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.1kHz의 균일 레이트로 변환합니다.

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 함수에서 사용하는 디폴트 방법입니다.

누락된 데이터가 있거나 입력값에 큰 간격이 있는 경우, 다른 보간 방법을 선택하여 누락된 데이터의 일부를 복구할 수 있습니다.

잡음이 적고 대역폭이 작은 신호의 경우 원래 신호를 재구성할 때 스플라인을 사용하면 매우 효과적입니다. 리샘플링하는 동안 3차 스플라인을 사용하려면 다음과 같이 '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를 선택하여 중간 그리드의 간격을 제어하면 됩니다.

약 3Hz 레이트로 진동하는, 부족감쇠 2차 필터의 계단 응답을 검토합니다.

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

계단 응답은 필터가 진동할 때는 높은 레이트로, 필터가 진동하지 않을 때는 낮은 레이트로 샘플링됩니다.

이제 디폴트 설정을 사용하여 100Hz로 신호를 리샘플링합니다.

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개 샘플의 그리드 샘플 레이트는 링잉 주파수의 나이퀴스트 샘플 레이트 6Hz보다 낮습니다.

그리드 샘플 레이트를 높이려면 정수 파라미터 p와 q를 지정하면 됩니다. resample이 그리드 샘플 레이트를 Q*Fs/P로 조정하고 신호를 보간한 다음, 내부 샘플 레이트 변환기를 적용하여(인자 P로 업샘플링하고 인자 Q로 다운샘플링) 원하는 샘플 레이트 Fs를 복원합니다. rat를 사용하여 p와 q를 선택합니다.

진동 링잉이 3Hz이기 때문에 샘플 레이트를 나이퀴스트 레이트보다 조금 높은 7Hz로 하여 그리드를 지정합니다. 1Hz의 여유 주파수는 감쇠 지수 포락선으로 인한 추가 주파수 성분이 차지합니다.

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

안티에일리어싱 필터 지정하기

다음 예제에서는 비행기 엔진의 스로틀 설정(Throttle Setting)을 측정하는 디지타이저의 출력값을 볼 수 있습니다. 스로틀 설정은 약 100Hz의 명목 레이트에서 불균일하게 샘플링되었습니다. 이 신호를 10Hz 레이트로 균일하게 리샘플링해 보겠습니다.

다음은 원래 신호의 샘플입니다.

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에 낮은 차단 주파수를 갖는 필터를 제공하여, 리샘플링된 신호를 더 매끄럽게 만들 수 있습니다.

우선, 명목 샘플 레이트가 약 100Hz가 되도록 그리드 간격을 설정합니다.

nominalFs = 100;

다음으로, 원하는 레이트를 얻는 데 적절한 p와 q의 값을 결정해야 합니다. 명목 레이트는 100Hz이고 원하는 레이트는 10Hz이기 때문에 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이 신호 경계를 벗어난 지점의 신호를 0이라고 가정하기 때문입니다. 이러한 불연속이 미치는 영향을 줄이려면, 신호의 양 끝점을 잇는 직선을 빼고 리샘플링을 수행한 뒤 그 직선을 원래 함수에 다시 추가하십시오. 이렇게 하려면 첫 번째 샘플과 마지막 샘플 사이를 잇는 직선의 오프셋과 기울기를 계산하고 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)

이제 추세 제거된 신호의 끝점 모두가 0에 가깝습니다. 이로써 과도 현상이 줄어듭니다. 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™ 문서를 참조하십시오.