fft matlab, scaling amplitude problem

Question

Bob GH 2013-5-30

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/77472-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 个评论
显示 -1更早的评论隐藏 -1更早的评论

Bob GH 2013-5-30

在 MATLAB Online 中打开

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])

请先登录，再进行评论。

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

Answer 1

Wayne King 2013-5-30

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/77472-fft-matlab-scaling-amplitude-problem#answer_87126

在 MATLAB Online 中打开

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 个评论
显示 2更早的评论隐藏 2更早的评论

Bob GH 2013-5-30

thanks for your reply Wayne King. i wrote a relevant part of my code

Bob GH 2013-5-30

在 MATLAB Online 中打开

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])

请先登录，再进行评论。

Answer 2

Azzi Abdelmalek 2013-5-30

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/77472-fft-matlab-scaling-amplitude-problem#answer_87124

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 个评论
显示无隐藏无

Bob GH 2013-5-30

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.