Finding the frequency value of a signal

510 次查看（过去 30 天）
Ramo
Ramo 2014-10-25
评论： Image Analyst 2021-10-26
How can I find the frequency of this signal?
Thanks,
  2 个评论
Yasemin G
Yasemin G 2021-10-26
Hello Ramo,
Can you share your .m file with me please I am having the same problem.
Image Analyst
Image Analyst 2021-10-26
@Yasemin G try this:
% Initialization steps.
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
imtool close all; % Close all imtool figures if you have the Image Processing Toolbox.
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 18;
%==================================================================================
% Image Analyst's code below:
period = 0.57e5
originalFrequency = 1/period
x = linspace(0, 2e5, 2000);
signalAmplitude = 0.85;
perfectSineWave = signalAmplitude * sin(2 * pi * (x - 0.08e5) / period);
subplot(2,1,1);
plot(x, perfectSineWave, 'b-');
grid on;
noiseAmplitude = 0.05;
yNoisy = perfectSineWave + noiseAmplitude * (2 * rand(1, length(x)) - 1);
hold on;
darkGreen = [0, 0.5, 0];
plot(x, yNoisy, '-', 'Color', darkGreen)
yline(0, 'Color', 'b', 'LineWidth', 2)
xlabel('Time or x', 'FontSize',fontSize);
ylabel('Signal', 'FontSize',fontSize);
title('Original Signal', 'FontSize',fontSize);
%==================================================================================
% Harry's code below:
% Assume we capture 8192 samples at 1kHz sample rate
Nsamps = 8192;
fsamp = 1000;
Tsamp = 1/fsamp;
t = (0:Nsamps-1)*Tsamp;
% Choose FFT size and calculate spectrum
Nfft = 1024;
[Pxx,f] = pwelch(yNoisy, gausswin(Nfft), Nfft/2, Nfft,fsamp);
% Plot frequency spectrum
subplot(2,1,2);
plot(f, Pxx, 'b-', 'LineWidth', 2);
ylabel('PSD', 'FontSize',fontSize);
xlabel('Frequency (Hz)', 'FontSize',fontSize);
grid on;
% Get frequency estimate (spectral peak)
[~,loc] = max(Pxx);
FREQ_ESTIMATE = f(loc)
caption = sprintf('Frequency estimate = %f Hz', FREQ_ESTIMATE);
title(caption, 'FontSize',fontSize);
g = gcf;
g.WindowState = 'maximized'

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Harry
Harry 2014-10-26
编辑：Harry 2014-10-26
Whenever you're interested in frequency content of a signal, the Fast Fourier Transform is often an excellent tool to use (see help fft). More specifically, Matlab's PWELCH function will provide a Power Spectral Density estimate using Welch's method:
[Pxx,F] = pwelch(X,WINDOW,NOVERLAP,NFFT,Fs)
Here is an example of how to use it to estimate frequency:
close all; clear all; clc;
% Assume we capture 8192 samples at 1kHz sample rate
Nsamps = 8192;
fsamp = 1000;
Tsamp = 1/fsamp;
t = (0:Nsamps-1)*Tsamp;
% Assume the noisy signal is exactly 123Hz
fsig = 123;
signal = sin(2*pi*fsig*t);
noise = 1*randn(1,Nsamps);
x = signal + noise;
% Plot time-domain signal
subplot(2,1,1);
plot(t, x);
ylabel('Amplitude'); xlabel('Time (secs)');
axis tight;
title('Noisy Input Signal');
% Choose FFT size and calculate spectrum
Nfft = 1024;
[Pxx,f] = pwelch(x,gausswin(Nfft),Nfft/2,Nfft,fsamp);
% Plot frequency spectrum
subplot(2,1,2);
plot(f,Pxx);
ylabel('PSD'); xlabel('Frequency (Hz)');
grid on;
% Get frequency estimate (spectral peak)
[~,loc] = max(Pxx);
FREQ_ESTIMATE = f(loc)
title(['Frequency estimate = ',num2str(FREQ_ESTIMATE),' Hz']);
  5 个评论
显示 3更早的评论
Image Analyst
Image Analyst 2017-4-6
f is returned by pwelch() - see the documentation for that.
x is defined as "x = signal + noise;" so I don't know why it would say that, unless you tried to use it before you defined it.

请先登录，再进行评论。

更多回答（2 个）

Image Analyst
Image Analyst 2014-10-25
I'd smooth it a bit with a 3rd order Savitzky-Golay filter, sgolayfilt() in the Signal Processing Toolbox, then I'd use findpeaks to get the period and 1/period is the frequency. Attached is a Savitzky-Golay filter demo.
  4 个评论
显示 2更早的评论
Ramo
Ramo 2014-10-26
>>
load uapp.txt
>> plot(uapp(:,1),uapp(:,2),'m-')
>> x = uapp(:,2);
>> smtlb = sgolayfilt(x,3,41);
subplot(2,1,1)
plot(1:200000, x(1:200000))
axis([0 200000 -15000 15000])
title('normal graph')
% grid
subplot(2,1,2)
plot(1:200000,smtlb(1:200000))
axis([0 200000 -15000 15000])
title('smoothed graph')
% grid
Image Analyst
Image Analyst 2020-3-27
Ramo, You chose a frame length, 41, that was far too small. It should be hundreds or thousands because you have 200 thousand data points (which is far too many to see on a single screen without zooming, by the way). Here is corrected code:
% Read in data.
data = dlmread('uapp.txt');
x = data(:, 1);
y = data(:, 2);
% Plot original noisy data.
subplot(1, 3, 2);
plot(x, y, '-');
title('Noisy Data', 'FontSize', 15);
grid on;
% Smooth the data
smoothedY = sgolayfilt(y, 4, 2001);
% Plot smoothed data.
subplot(1, 3, 3);
plot(x, smoothedY, 'r-', 'LineWidth', 2);
grid on;
title('Smoothed Data', 'FontSize', 15);
% Plot both original and smoothed data.
subplot(1, 3, 1);
plot(x, y, '-');
title('Noisy Data', 'FontSize', 15);
grid on;
hold on;
plot(x, smoothedY, 'r-', 'LineWidth', 2);
title('Both Noisy and Smoothed Data', 'FontSize', 15);
% Maximize the window.
g = gcf;
g.WindowState = 'maximized'

请先登录，再进行评论。

Ramo
Ramo 2014-10-26
I got 20000 points (i.e. x and y) which i load from a txt file. then I plot them using plot(x,y).. what is the next step? Thanks
  7 个评论
显示 5更早的评论
Ramo
Ramo 2014-10-30
If you could explain the equation you have used please? and by the way I have tested your method with another signal it seems that it gave me the right frequency, however the sine wave didnt match at all.
Rodrigo Picos
Rodrigo Picos 2020-3-27
Many years later.....
You should use 'sin3' or 'sin4', and check if you need to use the third or fourth component.
Hint: call F=fit( ) without the ending ;

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Spectral Estimation 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by