Main Content

이 번역 페이지는 최신 내용을 담고 있지 않습니다. 최신 내용을 영문으로 보려면 여기를 클릭하십시오.

신호 유사성 측정하기

이 예제에서는 신호 유사성을 측정하는 방법을 보여줍니다. 이는 다음과 같은 문제의 해결에 도움이 됩니다. 길이가 서로 다르거나 샘플링 레이트가 서로 다른 신호를 비교하려면 어떻게 해야 합니까? 측정값에 신호가 있는지 또는 잡음만 있는지 확인하려면 어떻게 해야 합니까? 두 신호가 관련이 있습니까? 두 신호 간의 지연을 측정하는 방법은 무엇입니까? 또한 이를 정렬하려면 어떻게 해야 합니까? 두 신호의 주파수 성분을 비교하려면 어떻게 해야 합니까? 유사성은 하나의 신호의 여러 부분에서 검출될 수도 있는 것으로, 그 신호가 주기성을 갖는지 여부를 결정하기도 합니다.

샘플링 레이트가 서로 다른 신호 비교하기

오디오 신호의 데이터베이스와 패턴 일치 응용 프로그램을 가지고 재생 중인 곡명을 알아내는 상황이라고 가정하겠습니다. 일반적으로 데이터는 메모리를 덜 차지하도록 낮은 샘플링 레이트로 저장됩니다.

load relatedsig.mat

figure
ax(1) = subplot(3,1,1);
plot((0:numel(T1)-1)/Fs1,T1,'k')
ylabel('Template 1')
grid on
ax(2) = subplot(3,1,2); 
plot((0:numel(T2)-1)/Fs2,T2,'r')
ylabel('Template 2')
grid on
ax(3) = subplot(3,1,3); 
plot((0:numel(S)-1)/Fs,S)
ylabel('Signal')
grid on
xlabel('Time (secs)')
linkaxes(ax(1:3),'x')
axis([0 1.61 -4 4])

Figure contains 3 axes. Axes 1 contains an object of type line. Axes 2 contains an object of type line. Axes 3 contains an object of type line.

첫 번째 서브플롯과 두 번째 서브플롯에서는 데이터베이스의 템플릿 신호를 보여줍니다. 세 번째 서브플롯은 데이터베이스에서 검색하려는 신호를 보여줍니다. 시계열만 봐도 신호가 두 템플릿 중 어떤 템플릿과도 일치하지 않습니다. 더 정밀하게 검사하면 신호가 실제로 길이와 샘플링 레이트에서 각각 다른 것을 알 수 있습니다.

[Fs1 Fs2 Fs]
ans = 1×3

        4096        4096        8192

길이가 다르므로 두 신호 간의 거리를 계산할 수 없지만, 이 문제는 신호의 공통 부분을 추출함으로써 손쉽게 해결할 수 있습니다. 사실, 길이를 항상 같게 만들 필요는 없습니다. 상호상관은 길이가 다른 신호 간에 수행할 수 있습니다. 단, 이 경우 샘플링 레이트는 반드시 동일해야 합니다. 이를 수행하기 위한 가장 안전한 방법은 더 낮은 샘플링 레이트로 신호를 리샘플링하는 것입니다. resample 함수는 리샘플링 과정 동안 안티에일리어싱(저역통과) FIR 필터를 신호에 적용합니다.

[P1,Q1] = rat(Fs/Fs1);          % Rational fraction approximation
[P2,Q2] = rat(Fs/Fs2);          % Rational fraction approximation
T1 = resample(T1,P1,Q1);        % Change sampling rate by rational factor
T2 = resample(T2,P2,Q2);        % Change sampling rate by rational factor

측정 데이터에서 신호 찾기

이제 xcorr 함수를 사용하여 신호 S와 템플릿 T1, T2의 상호상관을 분석하여 일치하는 부분이 있는지 확인할 수 있습니다.

[C1,lag1] = xcorr(T1,S);        
[C2,lag2] = xcorr(T2,S);        

figure
ax(1) = subplot(2,1,1); 
plot(lag1/Fs,C1,'k')
ylabel('Amplitude')
grid on
title('Cross-correlation between Template 1 and Signal')
ax(2) = subplot(2,1,2); 
plot(lag2/Fs,C2,'r')
ylabel('Amplitude') 
grid on
title('Cross-correlation between Template 2 and Signal')
xlabel('Time(secs)') 
axis(ax(1:2),[-1.5 1.5 -700 700 ])

Figure contains 2 axes. Axes 1 with title Cross-correlation between Template 1 and Signal contains an object of type line. Axes 2 with title Cross-correlation between Template 2 and Signal contains an object of type line.

첫 번째 서브플롯은 신호와 템플릿 1의 상관성이 낮음을 나타내는 반면, 두 번째 서브플롯의 높은 피크는 신호가 두 번째 템플릿에 존재함을 나타냅니다.

[~,I] = max(abs(C2));
SampleDiff = lag2(I)
SampleDiff = 499
timeDiff = SampleDiff/Fs
timeDiff = 0.0609

상호상관의 피크는 신호가 템플릿 T2에 존재하고 61ms 후 시작됨을 의미합니다. 다시 말해서, 신호 T2가 SampleDiff에 표시된 대로 499개 샘플만큼 신호 S보다 앞섭니다. 이 정보를 활용하여 신호를 정렬할 수 있습니다.

신호 간의 지연을 측정하고 정렬하기

여러 다른 센서에서 데이터를 수집하고, 다리 양쪽에서 자동차로 인해 발생하는 진동을 녹음하는 상황을 가정해 보겠습니다. 신호를 분석할 때 신호를 정렬해야 할 수 있습니다. 동일한 샘플링 레이트로 작동하는 3개의 센서가 있고 이 센서가 동일한 이벤트로 인해 발생하는 신호를 측정한다고 가정합니다.

figure
ax(1) = subplot(3,1,1);
plot(s1)
ylabel('s1')
grid on
ax(2) = subplot(3,1,2); 
plot(s2,'k')
ylabel('s2')
grid on
ax(3) = subplot(3,1,3); 
plot(s3,'r')
ylabel('s3')
grid on
xlabel('Samples')
linkaxes(ax,'xy')

Figure contains 3 axes. Axes 1 contains an object of type line. Axes 2 contains an object of type line. Axes 3 contains an object of type line.

finddelay 함수를 사용하여 두 신호 간의 지연을 검출할 수도 있습니다.

t21 = finddelay(s1,s2)
t21 = -350
t31 = finddelay(s1,s3)
t31 = 150

t21은 s2가 s1보다 350개 샘플만큼 지연됨을 나타내고, t31은 s3이 s1보다 150개 샘플만큼 앞서 있음을 나타냅니다. 이제, 이 정보를 활용하여 시간을 기준으로 신호를 이동하여 3개 신호를 정렬할 수 있습니다. 또한 alignsignals 함수를 직접 사용하여 신호를 정렬할 수 있습니다. 이 함수는 둘 중 더 빠른 신호를 지연시켜 두 신호를 정렬합니다.

s1 = alignsignals(s1,s3);
s2 = alignsignals(s2,s3);

figure
ax(1) = subplot(3,1,1);
plot(s1)
grid on 
title('s1')
axis tight
ax(2) = subplot(3,1,2);
plot(s2)
grid on 
title('s2')
axis tight
ax(3) = subplot(3,1,3); 
plot(s3)
grid on 
title('s3')
axis tight
linkaxes(ax,'xy')

Figure contains 3 axes. Axes 1 with title s1 contains an object of type line. Axes 2 with title s2 contains an object of type line. Axes 3 with title s3 contains an object of type line.

신호의 주파수 성분 비교하기

파워 스펙트럼은 각 주파수에 존재하는 전력을 표시합니다. 스펙트럼 일관성(Spectrum Coherence)은 신호 간의 주파수 영역 상관관계를 나타냅니다. 0 방향의 일관성 값은 대응하는 주파수 성분이 상관 관계가 없음을 나타내는 반면, 1 방향의 값은 대응하는 주파수 성분이 상관 관계가 있음을 나타냅니다. 두 신호와 각각 해당하는 파워 스펙트럼이 있다고 가정하겠습니다.

Fs = FsSig;         % Sampling Rate

[P1,f1] = periodogram(sig1,[],[],Fs,'power');
[P2,f2] = periodogram(sig2,[],[],Fs,'power');

figure
t = (0:numel(sig1)-1)/Fs;
subplot(2,2,1)
plot(t,sig1,'k')
ylabel('s1')
grid on
title('Time Series')
subplot(2,2,3)
plot(t,sig2)
ylabel('s2')
grid on
xlabel('Time (secs)')
subplot(2,2,2)
plot(f1,P1,'k')
ylabel('P1')
grid on
axis tight
title('Power Spectrum')
subplot(2,2,4)
plot(f2,P2)
ylabel('P2')
grid on
axis tight
xlabel('Frequency (Hz)')

Figure contains 4 axes. Axes 1 with title Time Series contains an object of type line. Axes 2 contains an object of type line. Axes 3 with title Power Spectrum contains an object of type line. Axes 4 contains an object of type line.

mscohere 함수는 두 신호 간의 스펙트럼 일관성을 계산합니다. 이 함수의 결과를 통해 sig1과 sig2가 약 35Hz와 165Hz의, 2개의 상관된 성분을 가진다는 것을 확인할 수 있습니다. 스펙트럼 일관성이 높은 주파수에서, 상관된 성분 간의 상대적인 위상은 교차 스펙트럼 위상으로 추정될 수 있습니다.

[Cxy,f] = mscohere(sig1,sig2,[],[],[],Fs);
Pxy = cpsd(sig1,sig2,[],[],[],Fs);
phase = -angle(Pxy)/pi*180;
[pks,locs] = findpeaks(Cxy,'MinPeakHeight',0.75);

figure
subplot(2,1,1)
plot(f,Cxy)
title('Coherence Estimate')
grid on
hgca = gca;
hgca.XTick = f(locs);
hgca.YTick = 0.75;
axis([0 200 0 1])
subplot(2,1,2)
plot(f,phase)
title('Cross-spectrum Phase (deg)')
grid on
hgca = gca;
hgca.XTick = f(locs); 
hgca.YTick = round(phase(locs));
xlabel('Frequency (Hz)')
axis([0 200 -180 180])

Figure contains 2 axes. Axes 1 with title Coherence Estimate contains an object of type line. Axes 2 with title Cross-spectrum Phase (deg) contains an object of type line.

35Hz 성분 간의 위상 지연은 -90도에 가까우며, 165Hz 성분 간의 위상 지연은 -60도에 가깝습니다.

신호에서 주기성 찾기

겨울 동안 사무실 건물에서 측정한 일련의 온도 값이 있다고 가정하겠습니다. 온도 값은 약 16.5주 동안 30분 간격으로 측정되었습니다.

load officetemp.mat  

Fs = 1/(60*30);                 % Sample rate is 1 sample every 30 minutes
days = (0:length(temp)-1)/(Fs*60*60*24); 

figure
plot(days,temp)
title('Temperature Data')
xlabel('Time (days)')
ylabel('Temperature (Fahrenheit)')
grid on

Figure contains an axes. The axes with title Temperature Data contains an object of type line.

화씨 70도대의 저온에서는 신호의 작은 변동을 분석하기 위해 평균을 제거해야 합니다. xcov 함수는 상호상관을 계산하기 전에 신호의 평균을 제거합니다. 이 함수는 교차공분산(Cross-covariance)을 반환합니다. 교차공분산을 정확히 추정하기 위해 최대 지연을 신호의 50%로 제한합니다.

maxlags = numel(temp)*0.5;
[xc,lag] = xcov(temp,maxlags);         

[~,df] = findpeaks(xc,'MinPeakDistance',5*2*24);
[~,mf] = findpeaks(xc);

figure
plot(lag/(2*24),xc,'k',...
     lag(df)/(2*24),xc(df),'kv','MarkerFaceColor','r')
grid on
xlim([-15 15])
xlabel('Time (days)')
title('Auto-covariance')

Figure contains an axes. The axes with title Auto-covariance contains 2 objects of type line.

자기공분산(Auto-covariance)에서 주요 변동과 작은 변동을 확인해 보십시오. 주요 피크와 작은 피크는 등거리로 나타납니다. 이를 확인하기 위해 피크 사이의 거리를 계산하고 플로팅합니다.

cycle1 = diff(df)/(2*24);
cycle2 = diff(mf)/(2*24);

subplot(2,1,1)
plot(cycle1)
ylabel('Days')
grid on
title('Dominant peak distance')
subplot(2,1,2)
plot(cycle2,'r')
ylabel('Days')
grid on
title('Minor peak distance')

Figure contains 2 axes. Axes 1 with title Dominant peak distance contains an object of type line. Axes 2 with title Minor peak distance contains an object of type line.

mean(cycle1)
ans = 7
mean(cycle2)
ans = 1

작은 피크(Minor Peak)는 주 7회, 주요 피크(Dominant Peak)는 주 1회 나타나고 있습니다. 이는 온도가 조절되는 건물에서 일주일 동안 측정된 데이터를 가져온 것이란 것을 생각해보면 이해할 수 있습니다. 첫 번째 7일 주기는 건물의 온도에 주 단위의 주기성이 있어 온도가 주말 동안 낮았다가 주중에 다시 정상으로 돌아오는 것을 나타냅니다. 1일 주기 동작은 일 단위의 주기성이 있어 야간에는 온도가 낮았다가 주간에는 온도가 높아지는 것을 나타냅니다.

참고 항목

| | | | | |