Experiment with combining signals to increase SNR.
I will try to add the signals together by computing their phase differences programmatically, shifting them, and adding them together.
First, some setup
foo = csvread('sig11_tag150.2_J20062.csv'); Fs = 1000000; % Sampling Frequency % Construct a 200KHz BPF N = 10; % Order Fpass1 = 198000; % First Passband Frequency Fpass2 = 202000; % Second Passband Frequency Apass = 1; % Passband Ripple (dB) h = fdesign.bandpass('N,Fp1,Fp2,Ap', N, Fpass1, Fpass2, Apass, Fs); Hd = design(h, 'cheby1'); samples = 670000:710000; % Filter the four channels out1 = filter(Hd, foo(samples,2)); out2 = filter(Hd, foo(samples,3)); out3 = filter(Hd, foo(samples,4)); out4 = filter(Hd, foo(samples,5));
Now compute the cross-correlation, based on the first signal
N = xcorr(out1, out2); x1 = fliplr(N(1:length(out1))'); N = xcorr(out1, out3); x2 = fliplr(N(1:length(out1))'); N = xcorr(out1, out4); x3 = fliplr(N(1:length(out1))');
Now work out the required shift
[a, shift2] = max(x1(1:30)); [a, shift3] = max(x2(1:30)); [a, shift4] = max(x3(1:30));
Then shift those signals
outs1 = out1 % First signal is not shifted, but give it a consistent name. outs2 = vertcat(out2(shift2:length(out2)), zeros(shift2-1,1)); outs3 = vertcat(out3(shift3:length(out3)), zeros(shift3-1,1)); outs4 = vertcat(out4(shift4:length(out4)), zeros(shift4-1,1)); origsig = [out1(15000:15100), out2(15000:15100), out3(15000:15100), out4(15000:15100)]; newsig = [out1(15000:15100), outs2(15000:15100), outs3(15000:15100), outs4(15000:15100)]; % Compare these two plots to check that they're now coherent figure subplot(2,1,1) plot(origsig) axis tight title('Original signals') subplot(2,1,2) plot(newsig) title('Shifted signals') axis tight % Now examine the addition of all of those newsig = out1+outs2+outs3+outs4;
Now we'll look at the SNR of the individual signals and the combined one
N = length(out1); k=-(N-1)/2:(N-1)/2; T=N/Fs; freq=k/T; figure subplot(2,1,1) bleh=fft(outs1); plot(freq, abs(bleh)); Pyy = bleh .* conj(bleh); Pyy = 10 * log10(Pyy); plot (freq, Pyy) subplot(2,1,2) bleh=fft(newsig); plot(freq, abs(bleh)); Pyy = bleh .* conj(bleh); Pyy = 10 * log10(Pyy); plot (freq, Pyy)
Well.. it looks like the signal has gone to a maximum of 22dB (note, not sure this scale makes sense), up from 0 (or from 15.1 dB in out3) and the noise is also about 10dB lower, so an SNR improvement of around 30 dB!
How about if the shift is applied to the unfiltered signal?
wide1 = foo(samples,2); wide2 = foo(samples,3); wide3 = foo(samples,4); wide4 = foo(samples,5); wide2 = vertcat(wide2(shift2:length(wide2)), zeros(shift2-1,1)); wide3 = vertcat(wide3(shift3:length(wide3)), zeros(shift3-1,1)); wide4 = vertcat(wide4(shift4:length(wide4)), zeros(shift4-1,1)); widesig = wide1 + wide2 + wide3 + wide4; figure subplot(2,1,1) bleh=fft(wide3); plot(freq, abs(bleh)); Pyy = bleh .* conj(bleh); Pyy = 10 * log10(Pyy); plot (freq, Pyy) title('Wide3') subplot(2,1,2) bleh=fft(widesig); plot(freq, abs(bleh)); Pyy = bleh .* conj(bleh); Pyy = 10 * log10(Pyy); plot (freq, Pyy) title('Combined wide')
This results in a signal increase from 9.719 dB (wide1) to 23.73 dB and noise increase from -16 dB to -9dB. The best signal received was 16.21dB on wide3, so this is still not bad (though I'm not sure it would ever make sense to combine unfiltered signals)
Out of interest, how does this approach compare to naively mashing uncorrected signals together?
dodgysig = foo(samples,2) + foo(samples,3) + foo(samples,4) + foo(samples,5); bleh=fft(outdodgy); plot(freq, abs(bleh)); Pyy = bleh .* conj(bleh); Pyy = 10 * log10(Pyy); plot (freq, Pyy) title('Unshifted combined signals')
Gives a maximum signal of 17.32 dB and noise of -9dB.
So just combining them is slightly better than just using the best signal alone, and 6dB worse than fiddling around with all of this. Or thereabouts. It might be interesting to combine them all and then filter them:
outdodgy = filter(Hd, dodgysig);
Max signal of 16.2 dB, and much less noise (< -60dB).
In summary:
Could next work on samples with multiple blips and try to pull them out by examining an fft output for maxima over a detection threshold.