您接收到一些数据,并想知道它是否与您测量的较长的流相匹配。即使数据被噪声破坏,相互关系也能让你做出这样的判断。
将一个环在桌面上旋转的记录加载到工作区中。截取一个一秒钟的片段,然后听一听。
load('Ring.mat') Time = 0:1/Fs:(length(y)-1)/Fs; m = min(y); M = max(y); Full_sig = double(y); timeA = 7; timeB = 8; snip = timeA*Fs:timeB*Fs; Fragment = Full_sig(snip); % To hear, type soundsc(Fragment,Fs)
绘制信号和片段。突出显示片段端点以供参考。
plot(Time,Full_sig,[timeA timeB;timeA timeB],[m m;M M],'r--') xlabel('Time (s)') ylabel('Clean') axis tight
plot(snip/Fs,Fragment) xlabel('Time (s)') ylabel('Clean') title('Fragment') axis tight
计算并绘制完整信号和片段的相互关系。
[xCorr,lags] = xcorr(Full_sig,Fragment); plot(lags/Fs,xCorr) grid xlabel('Lags (s)') ylabel('Clean') axis tight
互相关最大的滞后是信号起始点之间的时间延迟。重新绘制信号,覆盖片段。
[~,I] = max(abs(xCorr)); maxt = lags(I); Trial = NaN(size(Full_sig)); Trial(maxt+1:maxt+length(Fragment)) = Fragment; plot(Time,Full_sig,Time,Trial) xlabel('Time (s)') ylabel('Clean') axis tight
重复这一过程,但在信号和片段中分别添加噪声。声音无法从噪音中分辨出来。
NoiseAmp = 0.2*max(abs(Fragment)); Fragment = Fragment+NoiseAmp*randn(size(Fragment)); Full_sig = Full_sig+NoiseAmp*randn(size(Full_sig)); % To hear, type soundsc(Fragment,Fs) plot(Time,Full_sig,[timeA timeB;timeA timeB],[m m;M M],'r--') xlabel('Time (s)') ylabel('Noisy') axis tight
尽管噪声水平很高,该程序仍能找到丢失的片段。
[xCorr,lags] = xcorr(Full_sig,Fragment); plot(lags/Fs,xCorr) grid xlabel('Lags (s)') ylabel('Noisy') axis tight
[~,I] = max(abs(xCorr)); maxt = lags(I); Trial = NaN(size(Full_sig)); Trial(maxt+1:maxt+length(Fragment)) = Fragment; figure plot(Time,Full_sig,Time,Trial) xlabel('Time (s)') ylabel('Noisy') axis tight