清空筛选条件
清空筛选条件

FFT is shifting my artificial signal?

4 次查看（过去 30 天）
Zach Ruble
Zach Ruble 2023-9-27
评论： Zach Ruble 2023-9-27
I'm working on a function that performs an FFT on a signal that comes from a weather station (The primary function of which is to separate a pressure signal into its mean, wave, and perturbation elements. This is refered to as a Reynold's Triple Decomposition).
The function is complete, but I need to validate the function. So I created an artificial signal that is comprised of a mean component, 3 dominant oscillations, and finally a random variable that will represent our perturbations. The artificial signal is generated using this code:
%% Inputs
n = 18000; % Number of values in signal
m = 2; % Range of random element
f1 = 0.01; % Frequency of primary resonance (Hz)
f2 = 0.2;
f3 = 3;
i1 = 0.5; % Intensity of wave components
i2 = 0.5;
i3 = 0.5;
avgVal = 500; % Mean value for signal
%% Setting up signal
X = 1:1:n;
Y_raw = ones(1,n);
Y_mean = Y_raw + (avgVal-1);
Y_wave = ((Y_raw .* (i1.*sin(f1.*X))) + (Y_raw .* (i2.*sin(f2.*X))) + (Y_raw .* (i3.*sin(f3.*X))));
Y_MpW = Y_wave + Y_mean;
Y_rand = -m + (m + m)*rand([1,n]);
Y_final = Y_MpW + Y_rand;
Then I have a second function that loads the data from the first as if it were a normal input.
%% Importing/Organizing data
load('Y_final.mat');
m_1 = length(Y_final);
%% Extraction of mean parameter (Step 1)
AvgElement = mean(Y_final);
DataMinusMean = Y_final - AvgElement;
%% Creation of the Hanning window
H_per = 0.05;
H_bound = m_1 * H_per;
H_NoTop = hann(H_bound*2);
H_Top = ones(1, m_1);
H_Top(1, 1:H_bound) = H_NoTop(1:H_bound);
H_Top(1, (m_1-H_bound):m_1) = H_NoTop(H_bound:end);
%% Period of the primary wave (Step 2)
Fs = 20; % Sampling frequency
nf = Fs/2; % Nyquist Frequency
T = 1/Fs; % Sampling period
L = m_1; % Length of signal
t = (0:L-1)*T; % Time vector
f = Fs*(0:(L/2))/L;
FT_P_dcs = fft(H_Top.*DataMinusMean);
Fu_P_dcs = ((FT_P_dcs./length(DataMinusMean)).*conj(FT_P_dcs./length(DataMinusMean)));
% conversion to the energy spectrum
n = (Fs.*(1:nf))./L;
DeltaN = Fs/L;
Ea_dcs = 2*Fu_P_dcs;
Spp_dcs = Ea_dcs./DeltaN;
P2_dcs = abs(Spp_dcs/L);
P1_dcs = P2_dcs(1:L/2+1);
P1_dcs(2:end-1) = 2*P1_dcs(2:end-1);
%% PLotting stuff
figure(1);
semilogx(f,(P1_dcs))
hold on;
xline(0.01, 'r-', 'LineWidth', 2);
xline(0.2, 'r-', 'LineWidth', 2);
xline(3, 'r-', 'LineWidth', 2);
grid on;
title("DCS", 'FontSize', 20);
xlabel("f (Hz)", 'FontSize', 20)
ylabel("Intensity", 'FontSize', 20)
hold off;
Now, the problem....
I have artificially generated a signal that has 3 distinct sinusoidal functions. One at 0.01 Hz, one at 0.2 Hz, and one at 3 Hz. When I plot my FFT (Along with 3 vertical lines at those corresponding points) it is very clear that there is something wrong. The FFT is picking up on 3 very strong frequencies but they are all about 3 times the expected frequency.
0.01 -> 0.3222
0.2 -> 0.6367
3.0 -> 9.5489
Here is an image of the plot:
Any sugestions for how I could go about fixing this?

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

采纳的回答

Paul
Paul 2023-9-27
Ran in:
Hi Zach,
a) create the X vector with spacing based on the sampling frequency, Fs
b) create the wave vector with sinusoids computed from 2*pi*f1, etc. because f1 is in Hz, not rad/sec
%% Inputs
n = 18000; % Number of values in signal
m = 2; % Range of random element
f1 = 0.01; % Frequency of primary resonance (Hz)
f2 = 0.2;
f3 = 3;
i1 = 0.5; % Intensity of wave components
i2 = 0.5;
i3 = 0.5;
avgVal = 500; % Mean value for signal
%% Setting up signal
%X = 1:1:n;
Fs = 20;
X = (0:n-1)/Fs;
Y_raw = ones(1,n);
Y_mean = Y_raw + (avgVal-1);
%Y_wave = ((Y_raw .* (i1.*sin(f1.*X))) + (Y_raw .* (i2.*sin(f2.*X))) + (Y_raw .* (i3.*sin(f3.*X))));
Y_wave = ((Y_raw .* (i1.*sin(2*pi*f1.*X))) + (Y_raw .* (i2.*sin(2*pi*f2.*X))) + (Y_raw .* (i3.*sin(2*pi*f3.*X))));
Y_MpW = Y_wave + Y_mean;
Y_rand = -m + (m + m)*rand([1,n]);
Y_final = Y_MpW + Y_rand;
m_1 = length(Y_final);
%% Extraction of mean parameter (Step 1)
AvgElement = mean(Y_final);
DataMinusMean = Y_final - AvgElement;
%% Creation of the Hanning window
H_per = 0.05;
H_bound = m_1 * H_per;
H_NoTop = hann(H_bound*2);
H_Top = ones(1, m_1);
H_Top(1, 1:H_bound) = H_NoTop(1:H_bound);
H_Top(1, (m_1-H_bound):m_1) = H_NoTop(H_bound:end);
%% Period of the primary wave (Step 2)
Fs = 20; % Sampling frequency
nf = Fs/2; % Nyquist Frequency
T = 1/Fs; % Sampling period
L = m_1; % Length of signal
t = (0:L-1)*T; % Time vector
f = Fs*(0:(L/2))/L;
FT_P_dcs = fft(H_Top.*DataMinusMean);
Fu_P_dcs = ((FT_P_dcs./length(DataMinusMean)).*conj(FT_P_dcs./length(DataMinusMean)));
% conversion to the energy spectrum
n = (Fs.*(1:nf))./L;
DeltaN = Fs/L;
Ea_dcs = 2*Fu_P_dcs;
Spp_dcs = Ea_dcs./DeltaN;
P2_dcs = abs(Spp_dcs/L);
P1_dcs = P2_dcs(1:L/2+1);
P1_dcs(2:end-1) = 2*P1_dcs(2:end-1);
%% PLotting stuff
figure(1);
semilogx(f,(P1_dcs))
hold on;
%xline(0.01, 'r-', 'LineWidth', 2);
%xline(0.2, 'r-', 'LineWidth', 2);
%xline(3, 'r-', 'LineWidth', 2);
grid on;
title("DCS", 'FontSize', 20);
xlabel("f (Hz)", 'FontSize', 20)
ylabel("Intensity", 'FontSize', 20)
hold off;

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Discrete Fourier and Cosine Transforms 的更多信息

标签

产品


版本

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by