How to transform time series into frequency dimain and obtain amplitude and phase

Question

AS 2020-7-14

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/564998-how-to-transform-time-series-into-frequency-dimain-and-obtain-amplitude-and-phase

评论： Rumana Afrin 2022-6-15

 clc;
clear all;
load data.dat;
time=0:2:3600;
time=time';
nn_ss=data(1:1800,1:88);
 for jj=1:88
    ss_nw = nn_ss(:,jj);
    Ts = mean(diff(time)); % sampling time
    Fs=500; %sampling frequency    ** It should be 1/Ts
    Fn=Fs/2; %Nyquist frequency
    L=length(time);
    Y=fft(ss_nw);
    fts=Y/L; %Normalized fft
    P_amp2=abs(Y/L);
    P_amp(:,jj) = P_amp2(1:L/2+1);
    P_amp(2:end-1,jj) = 2*P_amp(2:end-1,jj);
    f=Fs*(0:(L/2))/L;
    phase(:,jj)=angle(Y).*180/pi;
    phase1(:,jj)=phase(1:L/2+1);
    %phase1(:,jj)=phase(1:901,jj);
    ly=length(Y);
    fs=(-ly/2:ly/2-1)/ly*Fs;
 end

I have use this code for FFT of the time series with dimension 1800 by 88. But, when I have transformed it into frequency domain I got frequency domain half of the time series and its corresponding amplitude and phase values. But, it should same length to time series. But, according to code it is not coming. Please help me to correct the code

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

dpb 2020-7-14

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/564998-how-to-transform-time-series-into-frequency-dimain-and-obtain-amplitude-and-phase#answer_465779

在 MATLAB Online 中打开

There's really nothing wrong with the code; you have followed the example from the documentation for FFT and the following lines

    P_amp2=abs(Y/L);
    P_amp(:,jj) = P_amp2(1:L/2+1);
    P_amp(2:end-1,jj) = 2*P_amp(2:end-1,jj);
    f=Fs*(0:(L/2))/L;

return the positive frequency half of the estimated PSD and normalized it to match total power in the double-sided PSD (P_amp2) by the factor of 2 applied to all frequencies except DC and Fmax which are not replicated in the two-sided spectrum.

The frequency vector f is for the one-sided PSD estimate.

If you really want the two-sided PSD, then don't take the positive half elements in P_amp; just use P_amp2 and the corresponding double-sided frequency vector. At that point, for convenience since apparently wouldn't need them, could rename P_amp2 to P_amp and fs to f.

BTW, the FFT and other MATLAB routines are vectorized, you can compute all the results using the original array and dispense with the for loop...

nn_ss=data(1:1800,1:88);           % if size(data)==>1800,88, then this is not needed; use the data array
Ts = mean(diff(time));             % sampling time (by column)
Fs=500;                            %sampling frequency    ** It should be 1/Ts
Fn=Fs/2;                           %Nyquist frequency
L=length(time);
Y=fft(nn_ss);                      % FFT() by column
fts=Y/L;
P_amp=abs(Y/L);                    % 2-sided PSD estimate
phase=angle(Y).*180/pi;
ly=length(Y);
fs=(-ly/2:ly/2-1)/ly*Fs;