Fourier Transform - position to velocity derivative

Hi everybody,

I am trying to do position to velocity differentiation of discrete points along X-axis and carry out the frequency response of the signal by FFT.

Can someone please help me understand how to take the frequency axis in final plot?.. Can we do this by using ODE45 or by writing a function??.. Thank you

x=9000; % the number of data points along a-axis(position based)

t = 0:0.0015:300; %(0.0015 = sampling cycle in seconds)

dx = diff(x);

dt = diff(t);

dxdt = dx./dt;

fft_dxdt = fft(dxdt); % DFT 

fft_dxdt_magnitude = abs(fft_dxdt);

figure()

plot(fft_dxdt_magnitude)

Answer 1

[Mona Mahboob Kanafi]

Hello,

Your problem when assigning frequencies to fft results is only related to the shift of frequencies. The frequencies which you have assigned are corresponded to the results of matlab fftshift (fft after being centered in zero). Therefore you have to apply fftshift to your fft results as below:

T = 0.0017; %Time increment(in second)
t = (0 : T : (16.15-T))'; % 9500 samples
fs = 1/T; %Sample rate in Hertz
N = length(t); % Number of samples
TotalTime = T * N;

df = fs/N; %Frequency increment
% df = 1/ TotalTime;
f = -fs/2: df: fs/2-df; % shifted frequencies corresponding to matlab fftshift

FFT_V_x = fft(V_x); %Fourier Transform
FFT_V_x_magnitude = abs(fftshift(FFT_V_x)); %Absolute Values 
plot(f ,FFT_V_x_magnitude)

If you don't get my explanation, I again suggest you to read this link for further explanation:

Plotting the DTFT using the output of fft

And about your doubt on sampling frequency, since you haven't resample your data, your sampling frequency and therefore assigned frequecies will not change. If you later have problem that the assigned frequencies have higher number of samples than your Vx and Ax vectors, try extrapolating your velocity and acceleration, or simply add zero at the beginning and end to make them fit in size (this does not change the dominant frequency components).

Hope these helps,

-Mona

[Sunny Math]

Thank you very much for your reply. No, I don't think its the same what I am looking for. Please go through the below clear explanation of my problem.

Answer 2

[Sunny Math]

I have samples of values on position along X and Y axis w.r.t time. I want to find Velocity and acceleration values from position values. For example: For X axis:

T = 0.0017sec; %Time increment(samples for second)

t = [0:0.0017:16.15]'; %Time for 9500 samples(Transpose of Row Vector since X is a Column Vector)

fs = 1/T; %Sample rate in Hertz

X = (9500x1 double); %Position values along X-axis

L = length(X); %Length of Samples along X

dx = diff(X); %Derivative of X along axis

dt = diff(t); %Derivative of t

V_x = dx./dt; %Velocity

dt2 = (dt(1:end-1)+dt(2:end))/2; %Second derivative of t

A_x = diff(V_x)./dt2; %Acceleration

Once I have values of Velocity and Acceleration, I have carried out the FFT over the values,

FFT_V_x = fft(V_x); %Fourier Transform

FFT_V_x_magnitude = abs(FFT_V_x); %Absolute Values

df = fs/L; %Frequency increment

f = -fs/2:df:fs/2-df; %Create the frequency axis

plot(f,FFT_V_x_magnitude)

I am not sure about the frequency axis for Velocity and acceleration since it involves differentiation for both. Kindly suggest me the steps if I am wrong.

[Mona Mahboob Kanafi]

Hello, I have updated my answer. Hope it helps you.

[Sunny Math]

Perfect.. Thank-you very much.

And what about the first derivatives of position vector along X axis to get corresponding velocity and acceleration?.. I have used diff command to carry out the first and second derivatives as you can see in my program.. Is the right way to carry on differentiation or I can do it by ODE45 command?..Is ODE45 applicable to only equations?

[Mona Mahboob Kanafi]

Using diff function for your problem is correct and common. As you mentioned, ODE is for ordinary differential EQUATIONs and your problem is just a simple differentiation. Since in matlab there is no direct differentiation, using diff is the option here.

[Sunny Math]

Thank you for your answer. Well, I have carried out the FFT and I'm not convinced with FFT diagram for the values of my program. I have attached matlab file and Velocity.jpg for your reference. Please could you tell where I am going wrong. I have cross checked several times but still I am not getting whats wrong..

[Mona Mahboob Kanafi]

Hello,

I couldn't download your code, but your velocity PSD looks quite normal. I think what you want to see is in logarithmic scale which you must use:

loglog(f,FFT_V_x_magnitude)

This must solve the problem, but if you still can't get what you want, you must start applying window function (bartlett, tukey,Welch,...) to your signal before applying fft like this, since your signal is not periodic for fft:

V_x = V_x .* bartlett(length(V_x))';

Also, I assumed that you know what that high peak in zero frequency means. It is just related to the mean value of your signal and you can remove signal mean value to be zero before applying fft, or otherwise you can just ignore that high peak. Hope these help you.

[Sunny Math]

Ok. I have gone through some articles on FFT on matlab community. Please find an attachment for the sample Excel sheet. The samples are with respect to position over time. Two columns refer along X and Y axis respectively . Kindly carry on the successive derivatives for velocity and acceleration (say for example along X-axis) and following FFT. Please give your valuable feedback on your valuation. Please attach plots. Thank you very much.