Euler's Formula and FFT

Question

Seweryn 2023-11-2

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2042081-euler-s-formula-and-fft

评论： Sulaymon Eshkabilov 2023-11-4

Hi guys, from the .csv file, I performed an FFT analysis in Simulink and obtained the following equations (pdf file). Take this equation (screen 1) for example, can I expand it to this form (screen 2)? I know, I can use e^ix = cos(x) + isin(x) but i don't know how to do that. Thanks in advance!

Best regards,

Seweryn

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Sulaymon Eshkabilov 2023-11-2

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2042081-euler-s-formula-and-fft#answer_1345716

编辑：Sulaymon Eshkabilov 2023-11-2

在 MATLAB Online 中打开

tek0005.csv

If understood correctly, you are trying to find a non-linear fit of the large flacturations using sine and cosine functions. If so, here is how it can be attained:

% Data import

D = readmatrix('tek0005.csv');

D(1:15, :) = []; % Data cleaning

f1 = 2; % found using fft

Model = @(a, t)(a(1)*sin(2*pi*f1*t)+a(2)*cos(2*pi*f1*t));

t = D(:,1);

y1 = D(:,2); % Channel 1

a0 = [1 1]; % Initially estimated guess values for coeff a

Coeffs1 = nlinfit(t, y1, Model, a0);

a = Coeffs1; % Found fit model coeffs

y_fit = (a(1)*sin(2*pi*f1*t)+a(2)*cos(2*pi*f1*t));

figure(1)

plot(t, y1, 'r')

hold on

plot(t, y_fit, 'k--', 'LineWidth',2)

legend('Data: CH1', 'Fit Model')

ylabel('Channel 1')

xlabel('Time')

hold off

% Similarly, you can do for the other channel data

f1 = 2; % found from fft

t = D(:,1);

y2 = D(:,3);

% NB: mean of y2 is NOT "0". Thus, it should be removed from y2 or added to the fit model!

a0 = [1 1]; % Initially estimated guess values for coeff a

Coeffs = nlinfit(t, y2, Model, a0)

Coeffs = 1×2

-0.0335 -0.0016

a = Coeffs;

y_fit = (a(1)*sin(2*pi*f1*t)+a(2)*cos(2*pi*f1*t))+mean(y2);

figure

plot(t, y2, 'b')

hold on

plot(t, y_fit, 'k--', 'LineWidth',2)

legend('Data: CH2', 'Fit Model')

ylabel('Channel 2')

xlabel('Time')

hold off

3 个评论
显示 1更早的评论隐藏 1更早的评论

Seweryn 2023-11-4

在 MATLAB Online 中打开

Hi @Sulaymon Eshkabilov, I have the following code, can I get the Fourier equations like in the "screen2" for maybe 50 harmonics with correct coefficients and get correct fitting in plots? In this code, I get Fourier equations for only one harmonic component and only sin without cos part. Thanks in advance!

data = readmatrix('tek0005.csv', 'Delimiter', ';', 'HeaderLines', 21);
time = data(:, 1);
voltage = data(:, 2);
current = data(:, 3);
fs = 1 / (time(2) - time(1));
N = length(time);
frequencies = (0:N-1) * fs / N;
voltage_fft = fft(voltage);
current_fft = fft(current);
a1_voltage = abs(voltage_fft(2)) / N;
a1_current = abs(current_fft(2)) / N;
total_harmonic_distortion_voltage = sqrt(sum(abs(voltage_fft(2:end)).^2)) / abs(voltage_fft(2));
total_harmonic_distortion_current = sqrt(sum(abs(current_fft(2:end)).^2)) / abs(current_fft(2));
f1 = frequencies(2);
equation_voltage = @(t) a1_voltage * sin(2 * pi * f1 * t) + mean(voltage);
equation_current = @(t) a1_current * sin(2 * pi * f1 * t) + mean(current);
fprintf('Fourier equation for voltage: %.4f * sin(2 * pi * %.4f * t) + %.4f\n', a1_voltage, f1, mean(voltage));
fprintf('Fourier equation for current: %.4f * sin(2 * pi * %.4f * t) + %.4f\n', a1_current, f1, mean(current));
fprintf('For voltage: a1 = %.4f, THD = %.4f\n', a1_voltage, total_harmonic_distortion_voltage);
fprintf('For current: a1 = %.4f, THD = %.4f\n', a1_current, total_harmonic_distortion_current);
t = linspace(min(time), max(time), 1000); % Create 1000 time points
voltage_fit = equation_voltage(t);
current_fit = equation_current(t);
figure;
subplot(2, 1, 1);
plot(time, voltage, 'b', t, voltage_fit, 'r');
title('Voltage over time');
xlabel('Time [s]');
ylabel('Voltage [V]');
legend('Data', 'Fitting');
subplot(2, 1, 2);
plot(time, current, 'b', t, current_fit, 'r');
title('Current over time');
xlabel('Time [s]');
ylabel('Current [A]');
legend('Data', 'Fitting');
phase_voltage = angle(voltage_fft);
phase_current = angle(current_fft);
phase1_voltage = phase_voltage(2);
phase1_current = phase_current(2);
fprintf('Phase angle for voltage: %.4f radians\n', phase1_voltage);
fprintf('Phase angle for current: %.4f radians\n', phase1_current);