MATLAB: Not getting the correct phase using fft, code posted

fftharmonic functionsphase angle

Hi everyone,
I am trying to extract the correct phase angle from a forced damped harmonic oscillator using matlab ode45() function to solve it and then using matlab fft to do fft analysis and try to extract the phase angle from the fft. Code:
tspan=[0,50];
init=[-2;0];
step=1;
% call the solver
[t1,y1]=ode45(@f1,tspan,init);
%Find FFT
m = length(y1) % Window length
Y1=fft(y1(:,1),m);
n = fix(m);
plot(t1(1:n,1),y1(1:n,1),t1(1:n,1),ifft(Y1),'*')
%Find the phase angle
FT_power2 = abs(Y1).^2;
FT_phase2= (((angle(Y1))) * 180/pi);
[c2,i2] = max(FT_power2);
phase = FT_phase2(i2)
Problem is I do not think I am getting the correct angle (or I am completely not understanding what I am suppose to get). The function in question being solved by ode45 is:
yprime = -c*y(2)-k*y(1)+4*sin(2*t); where c,k,y(1) and y(2) are known.
Shouldn't I get angle of 45 (since it would 0 + the pi/2 lag from the sin function)? I am getting:
178.7608 in degrees.
Please help me out if possible I am at complete loss.
PS:I am completely new to fft and such.

Best Answer

You could get rid of the transient response by subtracting off the output of
yprime = -c*y(2)-k*y(1)  %no sinusoid
This should leave you with just the sinusoidal part of the response only. You should then be able to get its phase by looking at the output sinusoid's values at t=0, pi/2,pi,3*pi/2, etc...
Related Question