Clear Filters
Clear Filters

Trying to ubderstand the power distribution in fft plot

3 views (last 30 days)
clear all
close all
clc
L=10;
n=1.45;
c=2.9979e8;
dt = 6e-12;
T=10*2*L*n/c;
eps0=8.854e-12;
A=80e-12;
t = (-T/2/dt:1:T/2/dt)*dt;
Nt=round(T/dt);
fsine = 1e9;
vsine = 1;
phi = vsine*sin(2*pi*fsine*t);
EL1t=1.274e7*exp(1i*phi);
FP=fft(phi);
fs=1/dt/Nt;
Fs=(-1/dt/2:fs:1/dt/2-1);
figure
Z=plot(Fs,fftshift(abs(fft(EL1t/Nt).^2*2*n*c*eps0*A)));
As you see from the obtained fft plot , the peak of the graph at 0Hz is around 60W , but I am struggling to understand how the power is distributed throughout the plot.
The given input power is 100W I think. Shouldn't the central frequency at 0Hz be around 100W.
Where did the rest of power is what I am not understanding...
  2 Comments
David Goodmanson
David Goodmanson on 29 Apr 2024
Hello Yogesh,
the signal you are taking the fft of is proportional to
exp(i*const*sin(const*t))
is that your intent?

Sign in to comment.

Accepted Answer

Mathieu NOE
Mathieu NOE on 30 Apr 2024
hello
the total power of the fft spectrum equals ~100W
pow_spectrum_total = 99.9502
but your fft spectrum spread it on multiple fft bins , and you are focusing only on the peak , not on the sum
L=10;
n=1.45;
c=2.9979e8;
dt = 6e-12;
T=10*2*L*n/c;
eps0=8.854e-12;
A=80e-12;
t = (-T/2/dt:1:T/2/dt)*dt;
Nt=round(T/dt);
fsine = 1e9;
vsine = 1;
phi = vsine*sin(2*pi*fsine*t);
EL1t=1.274e7*exp(1i*phi);
FP=fft(phi);
fs=1/dt/Nt;
Fs=(-1/dt/2:fs:1/dt/2-1);
pow_spectrum = fftshift(abs(fft(EL1t/Nt).^2*2*n*c*eps0*A));
pow_spectrum_total = sum(pow_spectrum)
pow_spectrum_total = 99.9502
figure
Z=plot(Fs,pow_spectrum);

More Answers (1)

David Goodmanson
David Goodmanson on 30 Apr 2024
Edited: David Goodmanson on 30 Apr 2024
Hello Yogesh,
The fact that all the abs(peaks)^2 add up to the expected total is just Parseval's theorem and has really nothing to do with the cause of the peaks. There are a couple of things going on here. The first is that the chosen frequency fsine = 1e9 does not have an integral number of cycles within the time window, and so it is not periodic in the time window.
t = 6e-12;
T=10*2*L*n/c;
t = (-T/2/dt:1:T/2/dt)*dt;
Nt =round(T/dt); % Nt = 161224
dt*Nt*fsine % cycles
ans = 967.3440
Just that situation all by itself means that a single peak in the frequency domain is not possible and there will be side peaks that mess up the result. By this I mean if you do an fft on
exp(2*pi*i*fsine*t) (1)
you will get a main peak and several smaller side peaks. This is because fsine = 1e9 is not commensurate with the frequencies in the fft. The code below replaces the time and frequency arrays with some similar ones for which 1e9 does fit exactly, having 1000 cycles. Then for the fft of (1) you get a single peak at the correct ampitude and no unwanted side peaks at all. After that change the fft calculation goes on just as before.
The more fundamental aspect is that you are transforming a function of the form
exp(vsine*sin(2*pi*fsine*t))
and its fourier transform definitely does have sidebands. The following expression gives their magnitude
exp(i*z*sin(theta)) = J_0(z)
+ Sum{k=1 to inf} J_2k(z)*[exp(i*2*k*theta) + exp(-i*2*k*theta)]
+ Sum{k=1 to inf} J_2k+1(z)*[exp(i*(2*k+1)*theta) - exp(-i*(2*k+1)*theta)]
Using
z = vsine and theta = 2*pi*fsine*t
you can identify the harmonics, and after collecting the all constants into C you can get analytic predictions for the peaks, which are shown with the red circles in the plot.
L=10;
n=1.45;
c=2.9979e8;
eps0=8.854e-12;
A=80e-12;
N = 2e5;
delt = 5e-12;
t = (-Nt/2:Nt/2-1)*delt;
delf = 1/(Nt*delt);
f = (-Nt/2:Nt/2-1)*delf;
fsine = 1e9; % N*delt*fsine = 1000 cycles
vsine = 1;
phi = vsine*sin(2*pi*fsine*t);
EL1t=1.274e7*exp(1i*phi);
figure(1)
plot(f,fftshift(abs(fft(EL1t/N).^2*2*n*c*eps0*A)));
hold on
C = 1.274e7^2*2*n*c*eps0*A
P0 = C*besselj(0,vsine)^2
P1 = C*besselj(1,vsine)^2
P2 = C*besselj(2,vsine)^2
fpks = [-2 -1 0 1 2]*1e9;
plot(fpks,[P2 P1 P0 P1 P2],'or')
xlim(1e11*[-.1 .1])
grid on
hold off
******************** also, note that
P0 +2*P1 +2*P2
ans = 99.8724

Community Treasure Hunt

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

Start Hunting!