fft matlab, scaling amplitude problem

Hi I've faced a problem in my homework coding. please help me to solve it.
I wrote a code for a problem, which i have its results. but after taking fft function from input, the amplitude of output is twice more than expected result.
what are the reasons that might be caused this?
I really appreciate it.

1 commentaire

I wrote a part of my code related to fft :
clc; clear all;
fm=40; T=10; a=0; N=32768; %2^15 Cm=1.5; h=(T-a)/N; omega=2*pi*fm/T; hh=1/h; Iamp=2.25e3; V(1)=-65;
mhf(1)=mhfinf(V(1));
mhs(1)=mhsinf(V(1));
mNap(1)=mNapinf(V(1));
EL=leak(V(1),mhf(1),mhs(1),mNap(1));
t(1)=a;
Iinp(1)=a;
for j=1:N-1
Iinp(j+1)=Iamp*sin(omega*(t(j)^2));
%%%%%%%%% 4 step method of Runge-Kutta
K1=h*fftNapH(t(j),[V(j); mhf(j); mhs(j); mNap(j)],EL,Iinp(j+1)); deffv1=K1(1,1);mhf1=K1(2,1);mhs1=K1(3,1);mNap1=K1(4,1);
K2=h*fftNapH(t(j)+(h/2),[V(j)+(deffv1/2);mhf(j)+(mhf1/2);mhs(j)+(mhs1/2);mNap(j)+(mNap1/2)],EL,Iinp(j+1)); deffv2=K2(1,1);mhf2=K2(2,1);mhs2=K2(3,1);mNap2=K2(4,1);
K3=h*fftNapH(t(j)+(h/2),[V(j)+(deffv2/2);mhf(j)+(mhf2/2);mhs(j)+(mhs2/2);mNap(j)+(mNap2/2)],EL,Iinp(j+1)); deffv3=K3(1,1);mhf3=K3(2,1);mhs3=K3(3,1);mNap3=K3(4,1);
K4=h*fftNapH(t(j)+h,[V(j)+deffv3;mhf(j)+mhf3;mhs(j)+mhs3;mNap(j)+mNap3],EL,Iinp(j+1)); deffv4=K4(1,1);mhf4=K4(2,1);mhs4=K4(3,1);mNap4=K4(4,1);
V(j+1)=V(j)+1/6*(deffv1+(2*deffv2)+(2*deffv3)+deffv4); mhf(j+1)=mhf(j)+1/6*(mhf1+2*mhf2+2*mhf3+mhf4); mhs(j+1)=mhs(j)+1/6*(mhs1+2*mhs2+2*mhs3+mhs4); mNap(j+1)=mNap(j)+1/6*(mNap1+2*mNap2+2*mNap3+mNap4);
t(j+1)=t(j)+h;
end
V2=fft(V1,N); V3=V2(1:N/2); amp=abs(V3); ampb=amp(2:length(amp)); % remove DC part f= (0:N/2-1)*hh/N;
figure(2),plot(f,ampb) xlabel('Freq') ylabel('Amplitude') axis([0 50 0 2000])

Connectez-vous pour commenter.

Réponses (2)

Wayne King
Wayne King le 30 Mai 2013
You should always show your code:
fs = 1000;
t = 0:1/fs:1-1/fs;
L = length(x);
xdft = fft(x)/L;
plot(abs(xdft))
Exactly as I expect two peaks with amplitude 0.5
Or
xdft = 2*fft(x)/L;
xdft = xdft(1:length(x)/2+1);
plot(abs(xdft))

4 commentaires

Bob GH
Bob GH le 30 Mai 2013
Thanks for your reply. i tried it, if i apply (/L) the result is absolutely different with expected result, but if i remove (/L) the result is same as mine.
Wayne King
Wayne King le 30 Mai 2013
Again, you're not showing your code (not sure why). If you are using the fft() for amplitude estimation, you want to divide by the length of the input signal.
Bob GH
Bob GH le 30 Mai 2013
thanks for your reply Wayne King. i wrote a relevant part of my code
I wrote a part of my code related to fft :
clc; clear all;
fm=40; T=10; a=0; N=32768; %2^15 Tfft=10; Cm=1.5; h=(T-a)/N; omega=2*pi*fm/T; hh=1/h; Iamp=2.25e3; V(1)=-65;
mhf(1)=mhfinf(V(1));
mhs(1)=mhsinf(V(1));
mNap(1)=mNapinf(V(1));
EL=leak(V(1),mhf(1),mhs(1),mNap(1));
t(1)=a;
Iinp(1)=a;
for j=1:N-1
Iinp(j+1)=Iamp*sin(omega*(t(j)^2));
%%%%%%%%% 4 step method of Runge-Kutta
K1=h*fftNapH(t(j),[V(j); mhf(j); mhs(j); mNap(j)],EL,Iinp(j+1)); deffv1=K1(1,1);mhf1=K1(2,1);mhs1=K1(3,1);mNap1=K1(4,1);
K2=h*fftNapH(t(j)+(h/2),[V(j)+(deffv1/2);mhf(j)+(mhf1/2);mhs(j)+(mhs1/2);mNap(j)+(mNap1/2)],EL,Iinp(j+1)); deffv2=K2(1,1);mhf2=K2(2,1);mhs2=K2(3,1);mNap2=K2(4,1);
K3=h*fftNapH(t(j)+(h/2),[V(j)+(deffv2/2);mhf(j)+(mhf2/2);mhs(j)+(mhs2/2);mNap(j)+(mNap2/2)],EL,Iinp(j+1)); deffv3=K3(1,1);mhf3=K3(2,1);mhs3=K3(3,1);mNap3=K3(4,1);
K4=h*fftNapH(t(j)+h,[V(j)+deffv3;mhf(j)+mhf3;mhs(j)+mhs3;mNap(j)+mNap3],EL,Iinp(j+1)); deffv4=K4(1,1);mhf4=K4(2,1);mhs4=K4(3,1);mNap4=K4(4,1);
V(j+1)=V(j)+1/6*(deffv1+(2*deffv2)+(2*deffv3)+deffv4); mhf(j+1)=mhf(j)+1/6*(mhf1+2*mhf2+2*mhf3+mhf4); mhs(j+1)=mhs(j)+1/6*(mhs1+2*mhs2+2*mhs3+mhs4); mNap(j+1)=mNap(j)+1/6*(mNap1+2*mNap2+2*mNap3+mNap4);
t(j+1)=t(j)+h;
end
V2=fft(V1,N); V3=V2(1:N/2); amp=abs(V3); ampb=amp(2:length(amp)); % remove DC part f= (0:N/2-1)*hh/N;
figure(2),plot(f,ampb) xlabel('Freq') ylabel('Amplitude') axis([0 50 0 2000])

Connectez-vous pour commenter.

Azzi Abdelmalek
Azzi Abdelmalek le 30 Mai 2013

0 votes

That means that, before calculating the fft, you've made some errors, which we can't find, because you have not posted the code.

2 commentaires

Bob GH
Bob GH le 30 Mai 2013
thanks for your reply. actually, It has 10++ line of code and i do not know which part do you need? is there any distinguished reason for my problem? because, my result is exactly twice of expected result.
Wayne King
Wayne King le 30 Mai 2013
Modifié(e) : Wayne King le 30 Mai 2013
Does 10++ mean 12 lines? If it is a reasonable number, please post it all. Did you look at the code I posted below?

Connectez-vous pour commenter.

Catégories

En savoir plus sur Fourier Analysis and Filtering dans Centre d'aide et File Exchange

Question posée :

le 30 Mai 2013

Community Treasure Hunt

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

Start Hunting!

Translated by