How i can determine sinusoid frequency ?

Hi; Pleas can anyone explain how i can determine sinusoid frequency if I have a step response where the data form is two vectors one for the acceleration and the other is time, when i transfer the plot from time domain to frequency domain i use fft; however i am confused how i can specify f0 when i read

http://blogs.mathworks.com/steve/2014/09/23/sinusoids-and-fft-frequency-bins/

I will appreciated for any explanation. Thanks in advance

Answer 1

[Star Strider]

It could be difficult to get a dominant peak with the signal you posted, since it could be a heterodyne of several frequencies. You will have to take the fft (link) of your signal and use your judgement.

It will be easier if you subtract the mean of the time-domain signal before you take the fft (link). This will eliminate the ‘0 Hz’ peak.

EDIT — 01:48 UCT 2017 04 02

The true resonant frequencies of your system are determined by the eigenvalues of your ‘A’ matrix. For your simulation, they appear to be the eigenvalues of:

Aqcar = [0 1 0 0;-kus/mus -c/mus k/mus c/mus;0 -1 0 1;0 c/ms -k/ms -c/ms];

[Muna Shehan]

Thanks for your reply. Dear Star Strider; the only reason I asked from the beginning about f0 to get a good understanding about fft which I suppose to use it to transfer the Acc. data vector from time domain to frequency domain. So your updated replay mean that I can get the sinusoid frequencies of the system response by eigenvalues(Aqcar). Kindly, can you check that based on the posted code is correct. (the transformation to frequency domain) domain)

[Star Strider]

My pleasure.

Your Fourier transform code appears to be correct. I cannot read the axis labels and tick labels on the image you posted, so I cannot comment on it.

[Muna Shehan]

Thanks for your reply. so the posted code (for frequency domain) is okay for a damped sine wave in time domain.

[Star Strider]

My pleasure.

I do not see any posted code other than the blog post. It is of course correct. The code you posted in response to Image Analyst’s Answer appears correct to me.

[Muna Shehan]

Thanks again . I mean the code that I posted to Image Analysis's Answer.

[Star Strider]

My pleasure.

Answer 2

[Image Analyst]

Subtract the mean and take the fft using fft() or pwelch(). Then scan it and look for the max magnitude. You didn't supply data so I can't do much other than speculate some code:

spectrum = fft(yourSignal - mean(yourSignal));
[maxValue, indexOfMaxValue] = max(abs(spectrum));

[Muna Shehan]

Thanks for your reply. I have posted the code. I have the dynamic response (Acc.) in time domain for the dynamic system ( quarter car suspension system). Now, I try to transfer this response (Acc.) to frequency domain. How I can do it correctly by using fft?. So far I used the fft example.

https://www.mathworks.com/help/matlab/ref/fft.html I appreciated for any help to transfer my signal from time domain to frequency domain.

    ms = 325;               % 1/4 sprung mass (kg)
    mus = 65;               % 1/4 unsprung mass (kg)
    kus = 232.5e3;          % tire stiffness (N/m)
    grav = 9.81;            % acceleration of gravity (m/s^2)
    v = 12.5;                 % vehicle velocity (m/s)
    dt = 0.001;             % simulation time step  
    % Road profile 
    sinbump = @(L,A) A*(0.5+0.5*(sin(linspace(-pi/2,pi*3/2,25*L))));
    B2 = sinbump(12, 0.10);
    road_x=linspace(0,99.99,600);
    road_z(1:300)=B2;
    road_z(301:600)=0;
    % Simulation time
    tmax = 5;
    x_time = 0:dt:tmax;
    dx = road_x(2) - road_x(1);             % spacial step for input data

    z0dot = [0 diff(road_z)/dt];            % road profile velocity  

    k=15024.9;
    c=1568.2;

        dt2 = dx/v;                             % time step for input data
        x = v*x_time;                 % time/space steps to record output
       % Construct linear state space model 
        Aqcar = [0 1 0 0;-kus/mus -c/mus k/mus c/mus;0 -1 0 1;0 c/ms -k/ms -c/ms];
        Bqcar = [-1 0 0 0]'; Cqcar = [1 0 0 0; 0 1 0 0; 0 0 1 0; 0 0 0 1; 0 -1 0 1]; 
        Dqcar = 0;
        qcar = ss(Aqcar,Bqcar,Cqcar,Dqcar);
        x0 = [0 0 0 0]';                        % initial state
        uu = interp1(road_x,z0dot,x);umf=1;    % prepare simulation input

        % Simulate quarter car model
        y = lsim(qcar,uu*umf,x_time,x0);                      
        Acc1 = [0 diff(y(:,4))'/dt2];       % sprung mass acceleration
    %Time domain plot
    figure(1);
    plot(x_time,Acc1);

    xlabel('Time (sec)','fontsize',14);
    ylabel('Acceleration (m/s^2)','fontsize',14);
    title({'Sprung mass response to a 0.1 m high step','Vehicle velocity = 45 km, Ks = 15024.9 N/m^{2}, Cs=1568.2 Ns/m'})

    %Frequency domain analysis
    Fs=1/dt;
    Lenght = size(x_time,2);
    NFFT = 2^nextpow2(Lenght(1,1)-1); % Next power of 2 from length of b
    Y1 = fft(Acc1,NFFT)/(Lenght(1,1)-1);%Y_FREQ = fft(y_time,N2)*T
    f = Fs/2*linspace(0,1,NFFT/2+1);

    yPower1 = 2*abs(Y1(1:NFFT/2+1));       %magnitude
    ydB1 = 20*log10(abs(Y1(1:NFFT/2+1))); %power

    spectrum = fft(Y1 - mean(Y1));
    [maxValue, indexOfMaxValue] = max(abs(spectrum));

    figure(2);
    plot(f,yPower1);
    set(gca, 'XLim',[0  30])

    xlabel('Frequency (Hz)','fontsize',14)
    ylabel('Acc.','fontsize',14)

    figure(3);
    plot(f,ydB1);

    xlabel('Frequency (Hz)','fontsize',14)
    ylabel('Power (dB)','fontsize',14)