Get spectrum. Fourier Transform

Hello! I need to get Fourier spectrum. I think, I am close to the goal, but have some misunderstandings. Variable 'signal' contains the array of input data (50006 points), step of measure is 0.1, so I have a signal in physical space (signal(time)) and attempt to achieve the same in Fourier space (FourierTransform(frequency)).

Here is my try:

%%discretisation in physical space
step_time=0.1;           
T=step_time*length(signal);   % the whole time of measurements
time_=0.1:step_time:T;        % according to Nyquist theorem

%%disret. in Fourier space
f_step=1/T;
F_duration = 1/step_time; 
f_frequency = 0:     f_step:      F_duration;
f_frequency(end) = [];     

%%get spectrum and get it normalized
Fourier_trans = fft(signal);
N_=length(Fourier_trans);
a=(Fourier_trans.*conj(Fourier_trans))/N_; % amplitude and normalization

but it goes wrong. What did I do incorrect..? Probably in last lines? Thank you in advance!

Answer 1

[Star Strider]

See the R2015a documentation for fft. Note the code between the first (top) 2 plot figures.

See if this works:

step_time=0.1;
sampling_frequency = 1/step_time;
nyquist_frequency = sampling_frequency/2;

%%get spectrum and get it normalized
Fourier_trans = fft(signal);
N_ = length(Fourier_trans);
a = abs(Fourier_trans)/N_; % amplitude and normalization
frequency_vector = linspace(0, 1, fix(N_/2)+1)*nyquist_frequency;
idx_vct = 1:length(frequency_vector);

figure(1)
plot(frequency_vector, a(idx_vct)*2)
grid

Note — This is UNTESTED CODE. It should work.

[Ivan Volodin]

looks better.. Can you explain why my code does not go? here

a=(Fourier_trans.*conj(Fourier_trans))/N_

I got normalized spectrum and here

f_frequency = 0: f_step: F_duration;

I got x-values. It should work too..

[Star Strider]

It is difficult for me to follow your code.

This assignment:

a=(Fourier_trans.*conj(Fourier_trans))/N_

calculates the square of the amplitude (power) not the amplitude. However if you want power rather than amplitude, your code is correct.

This assignment:

f_frequency = 0: f_step: F_duration;

does not consider that the fft is symmetrical and that the way the fft function calculates it, only the first half of the vector (from 0 Hz to the Nyquist frequency) is what you want. Your ‘f_frequency’ assignment assigns values to the frequencies beyond the Nyquist frequency that do not actually exist in the fft output.

The highest uniquely resolvable frequency in a sampled signal is the Nyquist frequency. There are no frequencies higher than that.

[Ivan Volodin]

Thanks a lot! Do you know, how get normalization in the case of

a=(Fourier_trans.*conj(Fourier_trans))

just dividing by N_ does not work correct

[Star Strider]

My pleasure.

This assignment:

a = (Fourier_trans.*conj(Fourier_trans))

calculates power (the square of amplitude), so you would have to divide by ‘N_^2’ to normalise it correctly.

[Ivan Volodin]

I have tried, but it does not lead to correct result. Please check this code:

step_t=0.01; 
first_step=0;
last_step=1;
t=first_step:step_t:last_step;         

%%generate signal
y=22*sin(2*pi*4*t)+15*sin(2*pi*42*t) ;

%%plotting graph
subplot(3,1,1), plot(t,y); 
%subplot(2,1,2), plot(f_freq_new,Y_new);

%%MatLab Guide
% http://www.mathworks.com/help/matlab/ref/fft.html 
fourier=fft(y );     
N_=length(fourier);
f_ = (1 / step_t)  *   (   0:    (N_/2)   )    /    N_; 

a=(fourier.*conj(fourier))/(N_^2); % mormalization
a = a(1:N_ /2+1);                                       
a(2:end-1) = 2*a(2:end-1);  
subplot(3,1,2), plot(f_,a); 
subplot(3,1,3),  loglog(f_,a);

See.. what is incorrect? I took dicretization from matlab's example and normalized it as you said