how to smooth or filter the signal like this?
2 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
How to smooth the signal like this?
I want to do fourier transform of the data attached. I only need the sharp peak in the data, I don't know how to
filter the noise in the data. what about Hanning window or something else?
Your help would be highly appreciated.
Best regards
load('data.mat')
plot(t,dip)
0 Kommentare
Antworten (4)
Bora Eryilmaz
am 19 Dez. 2022
Bearbeitet: Bora Eryilmaz
am 19 Dez. 2022
You can apply a low-pass filter to remove the noise in the signal before taking the Fourier Transform of the data:
load('data.mat')
plot(t,dip)
hold on
Fs = 1/(t(2)-t(1));
[b,a] = butter(4, 0.03/(Fs/2), 'low');
y = filter(b,a, real(dip));
plot(t,y)
By the way, your original data, dip, seems to have imaginary values in it. You would probably want to get rid of those.
0 Kommentare
Star Strider
am 19 Dez. 2022
Bearbeitet: Star Strider
am 19 Dez. 2022
I saved that code.
I did not save the non-code text that went with it. I am not certain that I can reproduce the image that you posted, because I do not understand what it is or how it was calculated. I plotted the Fourier transform three times here, once with the absolute amplitude, once with the amplitude in decibels, and once with the log amplitude.
Reposted —
LD = load(websave('data','https://www.mathworks.com/matlabcentral/answers/uploaded_files/1236142/data.mat'));
t = LD.t(:);
dip = LD.dip(:);
figure
plot(t, real(dip))
grid
xlabel('t')
ylabel('AMplitude')
title('Original Signal')
xlim([min(t) max(t)])
Ts = mean(diff(t)); % Sampling Interval
Fs = 1/Ts; % Sampling Frequency
Fn = Fs/2; % Nyquist Frequency
L = numel(t);
NFFT = 2^nextpow2(L); % For Efficiency
FTdip = fft(dip-mean(dip).*hamming(L),NFFT)/L; % Windowed Fourier Transform
Fv = linspace(0, 1, NFFT/2+1)*Fn; % Frequency Vector
Iv = 1:numel(Fv); % Index Vector
figure
plot(Fv, abs(FTdip(Iv))*2, 'DisplayName','Data')
grid
xlabel('Frequency')
ylabel('Magnitude (Absolute)')
title('One-Sided Fourier Transform')
% xlim([0 0.5])
[maxPeak,idx] = max(abs(FTdip(Iv))*2)
maxPeakFreq = Fv(idx)
passband = maxPeakFreq*[0.75 1.25]
xline(passband, '-r', 'DisplayName','Passband')
legend('Location','best')
figure
plot(Fv, mag2db(abs(FTdip(Iv))*2), 'DisplayName','Data')
grid
xlabel('Frequency')
ylabel('Magnitude (dB)')
title('One-Sided Fourier Transform')
% xlim([0 0.5])
xline(passband, '-r', 'DisplayName','Passband')
legend('Location','best')
figure
semilogy(Fv, abs(FTdip(Iv))*2, 'DisplayName','Data')
grid
xlabel('Frequency')
ylabel('Magnitude')
title('One-Sided Fourier Transform')
% xlim([0 0.5])
xline(passband, '-r', 'DisplayName','Passband')
legend('Location','best')
filtered_dip = bandpass(real(dip), passband, Fs, 'ImpulseResponse','iir');
figure
plot(t, filtered_dip)
grid
xlabel('t')
ylabel('Amplitude')
title('Filtered Signal')
xlim([min(t) max(t)])
This code selects the maximum frequency and defines a passband based on it. It then filters the signal using that passband. (The passband is also depicted in the frequency domain plot to illustrate what that filter passes.) I do not suggest a narrower passband, however shifting the existing passband or enlarging it (or both) are straightforward.
EDIT — (19 Dec 2022 at 20:36)
Looking at thid further, I am not certain how to get the result you want from this isgnal. The Savitzky-Golay filter (sgolayfilt) may be the way to go on it (another option being a FIR comb filter), however it is not obvious to me how to get the result depicted in your latest plot image.
% new_dip = lowpass(real(dip), 0.1, Fs, 'ImpulseResponse','iir')
sgfilt_dip = sgolayfilt(real(dip), 3, 127);
figure
plot(t, sgfilt_dip)
grid
xlabel('t')
ylabel('Amplitude')
title('Savitzky-Golay Filtered Signal')
xlim([min(t) max(t)])
FTdip = fft(sgfilt_dip-mean(sgfilt_dip).*hamming(L),NFFT)/L; % Windowed Fourier Transform
figure
semilogy(Fv, abs(FTdip(Iv))*2, 'DisplayName','Data')
grid
xlabel('Frequency')
ylabel('Magnitude')
title('One-Sided Fourier Transform Of Savitzky-Golay Filtered Signal')
% xlim([0 0.5])
xline(passband, '-r', 'DisplayName','Passband')
legend('Location','best')
Experiment with the sgolayfilt function. Also, it would help if you described in some detail the result you want, since that is not obvious to me.
.
0 Kommentare
Image Analyst
am 19 Dez. 2022
Not sure what you consider a sharp peak. Is it anything above 0.25 or below -0.25? If so, you can just zero out everything with a magnitude less than 0.25:
load('data.mat')
dip(abs(dip) < 0.25) = 0;
plot(t, dip)
Then do whatever you want, such as filtering in the spectral domain.
0 Kommentare
Siehe auch
Kategorien
Mehr zu Digital Filtering finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!